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range  variability  of  the  ocean  bottom.  There  are  two  approximations  that 
arise  in  applications  of  coupled  mode  theory.  The  adiabatic  approximation, 
which  becomes  valid  in  the  limit  of  slight  range  variability  of  the  medium, 
involves  the  neglect  of  the  mode-mode  energy  coupling  process.  Another 
approximation  which  is  inherent  in  the  theory  involves  a  boundary  condition 
approximation  which  must  be  made  in  applications  involving  range  variable 
boundaries.  The  adiabatic  approximation  is  examined  with  respect  to  its 
validity  as  a  function  of  sediment  type  and  bottom  range  variability 
(bottom  slope  and  radial  sound  speed  gradient).  The  raypath  properties  of 
the  adiabatic  approximation  are  also  examined  with  respect  to  the  multipath 
conversion  process  caused  by  a  sloping  bottom  and  are  shown  to  agree  with  a 
ray  theory  view  of  the  process. ^The  boundary  condition  approximation  which 
arises  in  applications  of  the  theory  to  problems  with  nonhorizontal  boundaries 
is  also  investigated.  It  is  shown  that  the  effect  of  this  approximation  is 
an  acoustic  field  that  does  not  have  the  proper  energy  flow  characteristics. 

A  modification  of  the  theory  that  is  first  order  in  the  bottom  slope  is 
derived  which  alleviates  this  problem.  Numerical  calculations  for  the  field 
in  a  wedge  shaped  waveguide  are  presented  to  illustrate  the  effects  of  the 
boundary  condition  approximation  and  to  show  how  the  modifications  to  the 
theory  produce  the  desired  results.  Throughout  the  dissertation  the  mathemati-| 
cal  formalism  is  presented  in  a  form  amenable  to  numerical  calculations.  Also 
numerical  methods  for  implementing  the  theory  are  described  and  discussed. 
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CHAPTER  I 


INTRODUCTION 

The  analytical  description  of  sound  propagation  in  a 
range  changing,  ocean  environment  is  a  topic  of  fundamental 
interest  to  acousticians  because  the  horizontal  stratification  of 
acoustic  media  is  a  condition  that  is  rarely  realized  in  the 
physical  world.  In  underwater  sound  propagation  there  is  a  genuine 
need  for  the  calculational  ability  to  realistically  model  a  range 
variable  ocean  environment.  A  wealth  of  acoustical  data  already 
exists  which  cannot  be  analyzed  properly  without  taking  into 
account  the  range  variability  of  the  acoustic  medium. 

The  types  of  methods  for  modeling  sound  propagation  in  a 
range  dependent  medium  are  the  same  as  for  a  horizontally  invariant 
medium  and  may  be  generally  classified  as  ray  theoretical  or  wave 
theoretical.  Ray  theory  is  an  analytic  method  that  associates  the 
propagation  of  energy  between  source  and  receiver  positions  with 
particular  raypath  trajectories.  Wave  theory  methods,  on  the  other 
hand,  are  based  on  the  solution  of  the  wave  equation  satisfied  by 
the  acoustic  field. 

The  main  advantage  of  a  ray  theory  approach  to  the 
description  of  sound  propagation  is  in  its  physically  intuitive 
notion  of  propagation  along  raypaths.  In  many  situations  a  ray 
theory  viewpoint  is  useful  for  developing  an  understanding  of  the 
physical  processes  and  phenomenology  of  underwater  sound 
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propagation.  There  are  several  disadvantages  associated  with  ray 

theory,  however.  The  raypath  description  of  sound  propagation  in 

1.1-1  3 

ray  theory  is  based  on  a  WKB  *  '  asymptotic  solution  of  the 
acoustic  wave  equation,  and  as  such,  is  valid  only  in  the  limit  of 
high  frequencies  and  slowly  varying  media.  Also  associated  with 
the  theory  are  infinities  connected  with  caustics^'^  and  shadow 
zones.  In  addition,  in  most  implementations  of  ray  theory  methods 
in  range  dependent  propagation  media  (see  Refs.  1.5  and  1.6)  one  is 
restricted  in  the  level  of  detail  with  which  the  ocean  bottom  may  be 
described. 

There  have  in  the  past  been  many  applications  of  ray 

theory  methods  to  sound  propagation  in  a  range  variable  environ- 
1.7  1.8 

ment.  Weston  "  ’  '  has  used  ray  theory  methods  to  investigate 
horizontal  refraction  and  guided  propagation  in  weakly  range  depen- 
dent  media.  Smith  '  has  employed  ray  theory  methods  to  examine  the 
averaged  transmission  loss  for  range  dependent  sound  channels,  and 
Harrison^ '  ^  has  used  raypath  invariants  discovered  by  Weston  to 
describe  horizontal  ray  curvature  ef fee  ts  in  ocean  basins  and  troughs, 
and  near  seamounts.  Modifications  to  asymptotic,  ray  theory  solu¬ 
tions  of  the  wave  equation  for  a  range  dependent  medium  have  also 
been  proposed  (see  Refs.  1.11-1.14).  These  modifications,  for  the 
most  part,  have  involved  slightly  different  asymptotic  expansions 
of  the  solution  to  the  acoustic  wave  equation  in  an  effort  to  lower 
the  frequency  range  of  applicability  and  to  impart  more  wave-like 
properties  to  the  resulting  solution. 
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When  one  considers  a  wave  theory  approach  to  the 
numerical  modeling  of  sound  propagation  in  a  range  dependent  medium 
there  are  three  ways  to  proceed.  One  way  is  a  full  three- 
dimensional  numerical  solution  of  the  wave  equation.  Another  is 
through  the  use  of  the  parabolic  equation  method.  The  final  method 
is  the  subject  matter  of  this  dissertation,  the  coupled  mode  theory 
approach.  A  wave  theory  description  of  sound  propagation  is  desir¬ 
able  because  it  is  an  exact  treatment  of  the  problem  in  formulation, 
valid  for  all  frequencies.  Most  applications  of  wave  theory,  how¬ 
ever,  are  to  low  frequency  sound  propagation  because  the  numerical 
sampling  of  the  medium,  which  is  wavelength  dependent,  becomes  so 
fine  at  higher  frequencies  that  a  numerical  solution  is  not 
practical. 

Of  the  previously  mentioned  wave  theoretical  methods, 
only  two  have  been  applied  to  acoustic  wave  propagation,  the  para¬ 
bolic  equation  method  and  coupled  mode  theory.  In  the  parabolic 
1.15-1.22 

equation  *  *  method,  which  has  seen  wide  use  in  the  underwater 

acoustics  modeling  community,  the  Helmholtz  equation  for  the 
acoustic  field  is  approximated  by  a  differential  equation  which  is 
parabolic  in  form.  The  nature  of  the  parabolic  approximation, 
however,  is  such  that  the  field  so  obtained  is  not  valid  at  ranges 
near  the  source  position  and  involves  the  neglect  of  backscattered 
energy.  The  parabolic  equation  method  is  capable  of  describing 
range  changes  in  the  water  column  but  suffers  from  some  serious 
disadvantages  concerning  the  description  of  the  bottom.  The 
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numerical  methods  employed  in  the  solution  of  the  parabolic 
differential  equation  are  such  that  the  bottom  is  characterized  by 
a  critical  grazing  angle  below  which  the  bottom  is  perfectly 
reflecting  and  above  which  it  is  perfectly  absorbing.  Needless  to 
say,  this  restriction  is  serious  if  one  is  interested  in  the 
influence  a  realistic  ocean  bottom  has  on  sound  propagation  in  the 
water  layer. 

The  subject  matter  of  this  dissertation  is  the  application 
of  coupled  mode  theory  to  the  description  of  sound  propagation 
in  a  range  changing,  ocean  environment.  Coupled  mode  theory  was 

proposed  for  use  in  underwater  acoustic  propagation  applications  by 

1  23  1  24 

Pierce  *  in  1965  and  again  by  Milder  ’  in  1969.  Coupled  mode 

theory  is  a  formalism  built  around  conventional  normal  mode  theory. 
As  with  other  wave  theory  methods  coupled  mode  theory  is  computa¬ 
tionally  feasible  only  at  low  frequencies.  One  advantage  to  a 
normal  mode  description  of  sound  propagation  is  that  the  ocean 
bottom  is  included  as  part  of  the  propagation  from  the  beginning 
with  its  acoustic  properties  capable  of  being  described  in  detail. 
The  detailed  description  of  the  bottom  afforded  by  coupled  mode 
theory  is  an  important  advantage  when  considering  low  frequency 
sound  propagation  in  an  ocean  waveguide  since  bottom  attenuation 
effects  are  reduced  at  low  frequencies;  with  the  reduced  attenua¬ 
tion,  acoustic  energy  refracting  through  the  bottom  can  greatly 
influence  the  acoustic  field  in  the  water  column  in  a  manner  that  is 
sensitive  to  the  geoacoustic  properties  of  the  subbottom.  Hence, 


the  coupled  mode  theory  formalism  is  expected  to  be  most  useful  in 
acoustic  propagation  problems  involving  significant  bottom 
interac  tion. 

The  coupled  mode  theory  formalism  is  not  restricted  to 

acoustics  applications.  Very  similar  types  of  theories  have  been 

1.25 

applied  to  electromagnetic  wave  propagation.  Wait  and  Spies 

used  coupled  mode  theory  to  examine  mode  conversion  in  atmospheric 

1  26 

VLF  radio  wave  propagation.  Wexler  has  also  used  the  coupled 

mode  formalism  in  describing  microwave  propagation  in  waveguides 

having  discontinuities  along  the  path  of  propagation.  In  a  book  by 
1.27 

Unger  *  the  coupled  mode  theory  has  been  used  to  describe  guided 

propagation  ir<  fiber  waveguides. 

In  underwater  acoustics  applications,  coupled  mode  theory 

has  been  applied  to  a  variety  of  range  dependent  propagation  media. 
1  28-1  30 

S.  T.  McDaniel  '  has  applied  coupled  mode  theory  to  the 

description  of  mode-mode  coupling  induced  by  a  randomly  rough, 

1  31  1.32 

water-sediment  interface.  Tappert  and  Dozier  *  ’  "  have  also 

used  the  theory  to  investigate  the  fluctuation  statistics  of  the 
acoustic  field  in  an  ocean  waveguide  having  a  water  sound  speed 
profile  that  varies  randomly  as  a  function  of  depth  and  range. 

Applications  of  coupled  mode  theory  to  underwater  acous¬ 
tics  problems  of  a  more  deterministic  nature  have  been  relatively 
few.  Most  applications  of  the  theory  have  involved  the  adiabatic 
approximation  which  is  discussed  in  Chapter  III.  Within  the  adia¬ 
batic  approximation,  which  is  expected  to  be  valid  for  weakly  range 
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dependent  media,  the  mode-mode  coupling  effects  are  ignored. 
1.33 

Graves  et  al.  '  applied  the  adiabatic  approximation  of  coupled 

mode  theory  to  an  isovelocity,  wedge-shaped  model  of  an  oceanic 

waveguide.  Nagl  et  al.*’^  and  Graves  et  al.**"^  considered  the 

adiabatic  approximation  applied  to  realistic  ocean  waveguides 

1  36 

having  range  variable  sound  speed  profiles.  Chwieroth  et  al. 
applied  the  full  formalism  of  coupled  mode  theory  to  the  problem  of 
an  acoustic  waveguide  having  a  range  variable,  parabolic  sound 
speed  profile.  Variations  of  coupled  mode  theory  where  the  acous¬ 
tic  medium  is  partitioned  into  range  segments  in  which  the  acoustic 

•  1  37 

properties  are  range  independent  have  been  used  by  Waits  '  and 

1  •  38 

MacPherson  to  describe  sound  propagation  in  a  waveguide  of 

range  dependent  dimensions. 

Until  now,  the  applications  of  coupled  mode  theory  to 
deterministic  problems  have  either  focused  on  range  variations  of 
the  properties  of  the  water  layer  or  have  employed  the  adiabatic 
approximation.  No  applications  of  the  full  coupled  mode  theory 
have  been  made  in  problems  involving  lateral  variations  in  the 
bottom  or  range  variable  boundaries.  The  material  in  this  disser¬ 
tation  focuses  on  these  types  of  problems  and  is  concerned  with  the 
application  of  the  theory  to  low  frequency,  acoustic  propagation 
over  bottom  structures  whose  depth  and  acoustic  properties  vary 
with  range.  Coupled  mode  theory  is  especially  well  suited  to  the 
description  of  acoustic  propagation  over  range  variable  ocean 
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bottoms  because  of  the  detailed  manner  in  which  the  bottom  may  be 
described. 

The  purpose  of  this  dissertation  is  threefold.  First, 
it  presents  a  general  mathematical  description  of  coupled  mode 
theory  in  a  form  that  is  convenient  for  numerical  calculation. 
Second,  it  examines  the  adiabatic  approximation  with  regard  to  its 
validity  when  applied  to  propagation  over  range  dependent  ocean 
bottoms.  Third,  it  exposes  an  inconsistency  in  the  theory  when 
applied  to  propagation  in  waveguides  having  nonhorizontal  bound¬ 
aries.  The  inconsistency  is  related  to  approximations  on  the 
physically  proper  boundary  conditions  to  be  satisfied  by  the 
acoustic  field  and  it  is  shown  that  the  inconsistency  can  be 
remedied  by  modification  of  the  theory. 

The  originality  of  the  work  presented  in  this  disserta¬ 
tion  lies  in  its  focus  on  the  ocean  bottom.  All  previous  applications 
of  the  theory  have  considered  range  variations  in  the  water  layer 
and  have  not  considered  sloping  boundaries  at  all.  Since  many 
important  propagation  situations  involve  propagation  over  sloping 
bottoms  and  bottoms  of  range  variable  composition,  it  is  clearly 
important  that  the  theory  and  its  various  approximations  be  placed 
on  a  firm  theoretical  foundation  for  use  in  problems  such  as  these. 

The  work  presented  in  this  dissertation  is  divided  into 
six  chapters  of  which  this  is  the  first.  In  Chapter  II  the  mathe¬ 
matical  formalism  of  coupled  mode  theory  as  originally  proposed  by 
1  23  1  24 

Pierce  *  and  Milder  *  is  presented  in  detail.  In  Chapter  II  all 


» 


quantities  occurring  in  the  theory  are  expressed  in  a  form  con¬ 
venient  for  numerical  computation.  Also  included  in  Chapter  II  is  a 
discussion  of  how  coupled  mode  theory  is  implemented  in  acoustic 
propagation  problems  of  practical  interest  along  with  a  review  of 
mode  theory  for  range  invariant  waveguides. 

In  Chapter  III  the  adiabatic  approximation  to  coupled 

mode  theory  is  considered.  The  approximation  and  its  implications 

are  discussed  and  the  energy  flow  characteristics  of  the  field  so 

obtained  are  examined.  Next  the  sensitivity  of  the  adiabatic 

approximation  to  the  range  dependence  of  the  bottom  is  considered* 

Lateral  changes  in  the  sediment  sound  speed  profile  and  water- 

sediment  interface  slope  are  examined  using  an  adiabacity  criterion 

1.24 


developed  by  Milder. 


Finally  the  adiabatic  approximation  is 


examined  with  respect  to  the  process  of  multipath  conversion.  Itj  is 

I 

concluded  that  the  geometric  (raypath)  properties  of  acoustic 
propagation  are  described  quite  well  in  the  adiabatic  approxi¬ 
mation. 

In  Chapter  IV  application  of  coupled  mode  theory  to 
propagation  over  range  variable  ocean  bottoms  is  considered. 
During  the  course  of  research  it  was  discovered  that  boundary 
condition  approximations  which  are  made  in  problems  involving  range 
variable  interfaces  cause  the  field  so  obtained  to  have  physically 
improper  energy  flow  characteristics.  These  improper  energy  charac¬ 
teristics  are  a  consequence  of  the  inconsistency  of  the  boundary 
condition  approximations.  This  problem  has  not  risen  until  now 


"r**^*-i 
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because  all  previous  applications  of  the  theory  have  been  to  wave¬ 
guides  with  horizontal  boundaries.  This  problem  is  exposed  and 
discussed  in  detail  in  Chapter  IV.  Also  considered  in  Chapter  IV  is 
a  correction  to  coupled  mode  theory  which  remedies  the  inconsistent 
boundary  condition  approximation  and  energy  flow  characteristics  of 
the  field. 

Chapter  V  presents  numerical  calculations  to  verify  the 
assertions  and  proposed  corrections  to  coupled  mode  theory  given  in 
Chapter  IV.  The  waveguide  model  used  in  the  numerical  calculations 
is  an  isovelocity  wedge  having  a  pressure  release  surface  and  a 
rigid  sloping  bottom.  The  calculations  presented  in  Chapter  V  are 
the  first  reported  results  of  an  application  of  coupled  mode  theory 
to  a  waveguide  with  nonhorizontal  boundaries,  excluding  the  adia¬ 
batic  approximation.  Calculations  of  the  rate  of  energy  transport 
and  various  other  quantities  associated  with  the  theory  are 
presented.  The  numerical  methods  employed  in  the  calculations  are 
also  described. 

Chapter  VI,  the  final  chapter  of  the  dissertation,  is 
devoted  to  a  summary  of  the  results  and  conclusions  of  Chapters  II- 
V.  Also  included  in  Chapter  VI  is  a  discussion  of  potentially 
useful  and  interesting  applications  of  coupled  mode  theory  along 
with  suggestions  for  future  work. 


CHAPTER  II 


FUNDAMENTALS  OF  COUPLED  MODE  THEORY 

The  purpose  of  this  chapter  is  to  describe  the 
mathematical  formalism  of  coupled  mode  theory  and  to  derive 
mathematical  expressions  which  will  be  used  throughout  this 
report.  In  the  first  part  of  this  chapter,  normal  mode  theory 
for  a  horizontally  stratified  propagation  medium  will  be  reviewed. 
Next,  range  dependence  of  the  medium  will  be  introduced  and  the 
mathematical  formalism  of  coupled  mode  theory  derived.  The  chapter 
will  conclude  with  a  brief  discussion  of  how  coupled  mode  theory  is 
implemented. 

In  the  derivations  performed  in  this  chapter,  a 
cylindrical  coordinate  system  will  be  employed.  The  symmetry  axis 
(z-axis)  is  perpendicular  to  the  ocean  surface  and  extends  downward 
in  the  positive  z  direction  toward  the  ocean  bottom.  Azimuthal 
symmetry  about  the  z-axis  is  also  assumed. 

A.  Review  of  Normal  Mode  Theory  for  a  Range  Invariant 

Medium 

In  a  horizontally  stratified  medium,  the  variations  in 
the  properties  of  the  waveguide  are  restricted  to  the  z-direction. 
The  types  of  variations  of  the  medium  encountered  in  the  ocean  are 
depth  changing  sound  speed  and  density.  These  changes  can  be  of  a 
continuous  or  discontinuous  nature.  Of  course  range  dependence  in 
the  ocean  is  always  present,  though  it  is  in  general  much  weaker 
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than  the  depth  dependence.  Therefore  the  assumption  of  horizontal 
stratification  of  the  ocean  medium  is  always  an  approximation  to 
the  actual  physical  situation,  and  depending  on  the  degree  of  range 
variability,  this  approximation  can  be  quite  good  or  quite  poor. 

The  time  independent  wave  equation  to  be  solved  is 

V"(p(r,  z)  +  k“(z)iKr,z)  =  0  .  (2.1) 

In  Eq.  (2.1),  ^  is  the  time  harmonic  acoustic  velocity  potential  for 

a  fluid  medium.  The  quantity  k(z)  is  equal  to  w/c(z)  where  io  is  2tt 

times  the  frequency  and  c(z)  is  the  sound  speed  in  the  fluid  medium 

as  a  function  of  depth.  Equation  (2.1)  is  not  valid  when  continuous 

depth  variations  of  density  within  a  layer  are  present.  Rutherford 
2  1 

and  Hawker  ’  have  shown  that  continuous  density  variations  of  a 
magnitude  likely  to  be  found  in  the  oceanic  waveguide  are  of 
negligible  importance.  Therefore,  throughout  this  dissertation, 
continuous  density  variations  within  layers  will  be  ignored. 
Discontinuous  density  variations  across  an  interface  separating  two 
media  of  differing  acoustic  properties  do  not  affect  the  form  of 
(2.1)  and  are  allowed. 

A  solution  of  Eq.  (2.1)  in  terms  of  acoustic  normal  modes 

may  be  obtained  in  a  number  of  ways  (see  Refs.  2. 2-2. 6).  The  method 

2.2 

employed  here  is  one  used  by  Mitchell.  '  Equation  (2.1)  can  be 
expressed  in  terms  of  products  of  radial  and  depth  functions  each 
depending  only  on  r  and  z,  respectively. 
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4»(r,z)  =  R(r)<j>(z) 


(2.2) 


The  differential  equation  governing  R(r)  and<Mz)  may  be  obtained 
by  substituting  Eq.  (2.2)  into  (2.1)  to  get 


~~2  +  y  +  ~2  +  k2(2)]  R(r)<Mz)  =  0 
idr  r  dT  dzZ  J 


(2.3) 


Equation  (2.3)  can  be  rewritten  as 


z)  .  (2.4) 


Since  the  left-hand  side  (LHS)  of  (2.4)  is  a  function  of  r  alone  and 

the  right-hand  side  (RHS)  of  (2.4)  a  function  of  z  alone,  the  LHS 

2 

and  RHS  must  both  be  equal  to  a  separation  constant  -<  .  Therefore 
the  equations  for  R( r)  and  <t>(z)  are  given  by 


—~2  R(r)  +  ~  —  R(r)  ♦  k2  R(r)  =  0 
dr 


(2.5) 


~2  <Kz)  +  |k2(z)  -  <2]  4>(z)  =  0 


(2.6) 


Equation  (2.5)  is  recognized  to  be  Bessel's  Equation.  If 
a  solution  involving  outgoing  waves  is  desired,  the  appropriate 
solution  to  (2.5)  for  e  lU)t  time  dependence  is 

R(r).>l<‘Vr)  ,2.7, 

where  ^  is  the  zero  order  Hankel  function  of  the  first  kind, 
o 

Equation  (2.6)  has  simple  analytic  solutions  only  for 

specific  types  of  depth  dependence  of  k(z).  For  example,  if  k(z)  is 

constant  in  depth,  <f)(  z)  is  a  combination  of  sine  and  cosine 

functions.  If  k  (z)  is  linear  with  depth  <J>(z)  is  given  by  combina- 

2.22.4  •  . 

tions  of  Airy  functions.  *  ’  "  The  boundary  conditions  to  be 
satisfied  by  (2.6)  are  obtained  from  the  boundary  conditions  to  be 
satisfied  by  the  total  field  ip  of  Eqs.  (2.1)  and  (2.2).  The 
boundary  conditions  on  and  its  derivatives  appropriate  for  an 
oceanic  waveguide  are  that:  (1)  <P  vanish  at  z=0  (pressure  release 
surface),  (2)  ^vanish  at  z=°°,  (3)  the  acoustic  pressure  be  continuous 
across  any  surfaces  of  discontinuity,  and  (4)  the  normal  component 
of  particle  velocity  be  continuous  across  any  surfaces  of  discon¬ 
tinuity.  Boundary  condition  (2)  is  a  consequence  of  the  fact  that 
the  sound  speed  structure  in  the  ocean  bottom  is  generally 
increasing  with  depth  in  the  bottom  and  attains  values  much  larger 
than  those  occurring  in  the  water  column.  Should  condition  (2)  not 
be  satisfied  guided  propagation  is  not  possible. 
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The  acoustic  pressure  p  is  given  by 


P  = 


P  |^  =  iwpip 


(2.8) 


Boundary  condition  (3)  requires  that  the  product  of  density  p  and 
depth  function  <Kz)  be  continuous  across  a  surface  of  discon¬ 
tinuity.  Boundary  condition  (4)  requires  that  the  normal  component 

3d) 

of  particle  velocity  given  by  be  continuous  across  any  surfaces 
of  discontinuity.  This  in  turn  requires  that 

9<t>(z) 

3z  (2.9) 


be  continuous  across  layer  boundaries. 

Figure  II. 1  shows  the  depth  structure  of  a  typical 
oceanic  waveguide.  The  first  layer  is  the  water  layer  which 
overlies  sediment  layers  which  are  treated  as  fluids  supporting 
only  compressional  waves.  The  water  and  sediment  layers  overlie  an 
infinite,  isovelocity  half-space  which  may  be  considered  to  have 
solid  properties  without  unduly  complicating  the  problem  (see 
Ref.  2.7).  The  material  densities  are  assumed  constant  within  each 
layer  while  the  compressional  wave  sound  speeds  are  allowed  to  vary 
continuously  in  the  water  and  sediment  layers. 
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In  terms  of  the  notation  of  Fig.  II.  1,  boundary  conditions 
(l)-(4)  become: 


(1) 

4,(0) 

=  0 

(2) 

<H°°) 

=  0 

(2.10) 

(3) 

Pi4>(l)(Hi) 

■  Pi.i 

(4) 

(2.11 ) 

The  boundary  conditions  (2.10)  and  (2.11)  along  with  Eq.  (2.6) 
define  an  eigenvalue  problem.  The  spectrum  of  eigenvalues  k^,  in 
general,  has  a  finite  number  of  discrete  elements  and  a  con¬ 
tinuously  distributed  component  (see  Refs.  2.8  and  2.9  and  Appendix 
A).  The  normal  mode  depth  equations  therefore  become 

-4  0n(z)  +  (k2(z)  -  kn2)  0n(z)  =  0  (2.12) 

dz 

subject  to  Eqs.  (2.10)  and  (2.11).  The  mode  functions  thus  deter¬ 
mined  form  an  orthogonal,  complete  set  of  functions  normalized  so 
that 

p<f>  (z)  <t>  (z)  dz  =  6  .  ( o  it-v 

n  m  n,m  iz.u; 

In  (2.13)  p  is  the  material  density  of  the  medium  taken  to  be 
constant  within  a  layer.  Any  type  of  compressi onal  wave,  acoustic 


17 


disturbance  in  the  oceanic  waveguide  can  be  described  in  terms  of 
the  radial  and  depth  functions  as  follows 

I 

iKr.z)  =  ]T  V^OV)  ^n(2)  *  (2.14) 

n 

The  £'  symbol  denotes  summation  over  discrete  and  integration  over 
continuum  eigenvalues.  In  most  applications  the  contribution  to 
Eq.  (2.14)  from  the  continuous  modes  is  negligible  (see  Ref.  2.2, 
pp.  24-25);  therefore,  as  is  customary  in  underwater  acoustics,  the 
continuum  contribution  to  the  field  will  be  neglected  in  the 
following  analysis.  See  Appendix  A  for  a  more  detailed  description 
of  the  properties  of  the  mode  depth  functions. 

One  particularly  useful  solution  of  the  wave  equation  is 
the  one  for  a  point  source  of  volume  velocity  located  at  ?=£  .  The 
equation  describing  the  field  in  this  situation  is 

V2iHr,z)  +  k2(z)tfi(r, z)  »  -4tt6(t-?o)  .  (2.15) 


After  Eq.  (2.15)  is  integrated  over  the  azimuthal  angle,  it  becomes 


/  92  13  32  \  , 


2  26(r-r  ) 

-H-)(Kr,z)  ♦  (z)tKr,z)  = - jr- 2  6(z-zq)(2.16) 


Upon  substituting  (2.14)  into  (2.16)  one  obtains 
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—  6(r-r  ) 6(z-z 
r  o 


where  Eq.  (2.12)  has  been  employed.  If  the  orthonormality 
ties  of  <f>  are  exploited  (see  Eq.  (2.13))  (2.17)  becomes 


,  dr" 


I  iL 

r  dr 


+  k 


H^(k  r)  =  —  6(r-r  )  p(z  )  <J>  (z  ) 
o  m  r  o  o  m  o 


The  constants  a^  can  be  determined  by  recognizing  that  Hq 
the  two-dimensional  freefield  Green's  function  (see  Ref 
Chapter  VII)  satisfying 


d  Id  .2 

— ~  +  —  -j—  +  k 
.2  r  dr  m 

dr 


t  6(r-r  ) 

H  1  (k  r)  =  - — 

o  m  7T  r 


Equations  (2.18)  and  (2.19)  combine  to  yield 


6(r-r  )  - 

— - - a  =  —  6(r-r  )  p(z  )  <p  (z  ) 

n  r  m  r  o  o  m  o 


a  =  itr  p(z  )  <p  (z  ) 
m  o  m  o 


Therefore  the  velocity  potential  due  to  a  point  source  at 
given  by 


>P(r,z) 


i7rp(zQ) 


(knr)  4>nCz)  *n( zQ) 


) (2.17) 

proper- 


.  (2.18) 

(k  r)  is 
m 

,  2.10, 

(2.19) 


(2.20) 

“k  • 

r=r  is 
o 


(2.21) 


» 
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where  the  contribution  from  the  continuum  has  been  neglected  (see 
Appendix  A).  With  the  derivation  of  Eq.  (2.21)  the  review  of  range 
stratified  normal  mode  theory  is  concluded.  The  next  section  of 
this  chapter  will  extend  the  ideas  of  normal  mode  theory  to  a  range 
dependent  medium. 


B.  Basics  of  Coupled  Mode  Theory 

The  theory  of  acoustic  propagation  via  coupled  normal 

.  .  .  2  1 1 

modes  was  proposed  independently  in  papers  written  by  Pierce  *  in 
2  12 

1965  and  Milder  '  in  1969.  The  mathematical  formalism  of  coupled 
mode  theory  was  described  in  these  two  papers  though  both  authors 
were  more  concerned  with  situations  in  which  the  mode  coupling 
process  could  be  ignored  (the  adiabatic  approximation). 

The  wave  equation  for  the  velocity  potential  in  a  medium 
with  depth  and  range  dependence  is 

V2if*(r , z)  +  k2(r,z)  iKr.z)  =  0  .  (2.22) 

In  Eq.  (2.22)  the  range  dependence  of  the  medium  enters  in  two  ways. 
One  is  through  the  range  dependence  of  the  wave  number.  Additional 
sources  of  range  dependence  enter  through  the  boundary  conditions. 
The  boundary  conditions  on  ip  are  applied  at  surfaces  of  discon¬ 
tinuity  whose  depth  might  vary  with  range  (i.e.,  sloping  layers). 
The  boundary  conditions  on  ijj  also  involve  the  density  which  could  be 
range  dependent.  The  range  dependence  of  the  boundary  conditions 
therefore  contributes  to  the  radial  dependence  of  the  field.  As  in 


~mr~ 


y 


20 


the  horizontally  stratified  case,  Eq.  (2.22)  is  expressed  in 

cylindrical  coordinates  and  azimuthal  symmetry  is  assumed.  The 

extention  to  three-dimensional  coupled  mode  theory  formulated  in 

cartesian  coordinates  is  given  later  in  this  chapter. 

Because  of  the  radial  dependence  ofk,  a  complete  separation 

2 

of  variables  is  possible  only  for  k  (r,z)  of  the  form 

k2(r,z)  =  k2(r)  ♦  k2(z) 

Since  this  form  for  the  wave  number  is  unduly  restrictive  and 
not  physically  realistic,  it  will  not  be  considered  further. 

Coupled  mode  theory  is  implemented  by  attempting  a 
partial  separation  of  variables.  In  analogy  with  conventional  mode 
theory,  Pierce  postulated  a  solution  of  Eq.  (2.22)  of  the  form 

>Kr,z)  =^Rn(r)  *n(z;r)  (2.23) 

n 

where  <J>n(z;r)  satisfies 

k2(r,z)  -  kn2(r)j<J>n(z;r)  =  0  .  (2.24) 

The  boundary  conditions  which  ^  is  to  satisfy  are  that  p<J>n  and 

3*  /9z  be  continuous  at  surfaces  of  discontinuity.  These 
n 

conditions  arise  from  requiring  continuity  of  pressure  and  z- 
component  of  particle  velocity  across  interfaces.  The  depth 
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functions  4>n  satisfy  (2.24)  locally  which  gives  rise  to  a  parametric 
dependence  of  the  mode  functions  and  eigenvalues  on  range. 

The  partial  separation  of  variables  indicated  by  (2.23) 
and  (2.24)  plus  boundary  conditions  is  exact  only  for  range 
variable  media  having  horizontal  interfaces  and  boundaries.  In 
this  case  the  z-component  of  particle  velocity,  8<)>/3z,  and  pressure 
must  be  continuous  across  boundaries.  Hence  the  boundary 
conditions  on  ^  may  be  transferred  directly  to  <t>n(z;r).  When 
nonhorizontal  boundaries  are  present,  it  is  the  normal  component  of 
particle  velocity,  V^<{)  ,  and  pressure  that  are  required  to  be 
continuous.  For  a  nonhorizontal  boundary  defined  by  z=H(r), 
is  given  by 


(2.25) 


h  -  . 

dr 

When  nonhorizontal  boundaries  are  present,  the  boundary  conditions 
on  \p  do  not  transfer  to  <f>n(2jr)  alone  and  nothing  is  gained  by 
expressing  the  field  as  in  Eq.  (2.23).  In  other  words,  not  even  a 
partial  separation  of  variables  is  possible,  greatly  complicating 
the  solution  of  Eq.  (2.22). 

The  intent  of  Pierce  was  that  the  formalism  of  coupled 
mode  theory  be  applied  to  problems  having  nonhorizontal  boundaries 
even  though  the  theory  is  not  strictly  valid  in  these  cases.  The 
reasoning  behind  this  assumption  was  that  for  small  slope  angles, 
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i.e.,  small  H,  the  normal  derivative  is  approx  imately  the  z- 
derivative;  hence,  the  true  boundary  conditions  could  be  replaced 
with  continuity  of  z  derivative  boundary  conditions  while  incurring 
negligible  errors  of  order  H.  It  turns  out  that  this  practice  is 
not  justified  and  can  introduce  serious  problems.  Chapter  IV  is 
concerned  with  correcting  coupled  mode  theory  for  use  in  problems 
with  nonhorizontal  boundaries.  The  rest  of  this  chapter  will 
derive  the  mathematical  formalism  of  coupled  mode  theory  as 
proposed  by  Pierce  and  will  not  be  concerned  with  the  improper 
boundary  conditions  that  must  be  imposed  in  problems  with  sloping 
boundaries. 


1 .  Derivation  of  Radial  Equations 

To  derive  the  differential  equation  for  the  radial 
functions  Rn(r)  Eq.  (2.23)  is  first  substituted  into  Eq.  (2.22)  to 
yie  Id 


R  <J>  R 

n  Tn  Yi 


*  ••  •  •  n  n  n 

>  R  4>  +  2R  <p  +  -2-^  +  —  <p  ♦  R  4> 

/  j  nyn  nYn  r  r  n  nrn 


(2.26) 


In  Eq.  (2.26)  and  in  the  rest  of  the  dissertation,  the  dot  symbol 
stands  for  derivative  with  respect  to  r.  The  depth  functions 
4>n(z,r),  at  each  range  point  r,  form  an  orthonormal  complete  set  of 
functions  such  that 
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3  4>mCz;r)  dz  =  6 
1,1  n  ,m 


(2.27) 


If  Eq.  (2.26)  is  multiplied  by  P <P  (z;r)  and  integrated  over  depth, 

m 

the  following  differential  equation  for  R  ( r)  results 

n 


R  (r)  +  —  K  (r)  + 
m  r  m 


k  2(r)  R  =  -  V  A  R  +  B  f—  +  2R  \(2  28) 
m  m  /  *  mn  n  mn  \  r  n/'**^ 


The  A  and  B  are  functions  of  range  defined  by 
mn  mn  J 


(r)  =  /  pU)  <l>m(z;r)  <t>n(z;r)  dz 

*o 


(2.29) 


B  (r)  =  f  p(z)  <}>  (z;r)  4>  (z;r)  dz 

mn  J  m  n 

o 


(2.30) 


and  are  referred  to  as  the  coupling  coefficients.  For  reasons 

related  to  numerical  solutions  of  differential  equations,  a  more 

convenient  differential  equation  is  sometimes  the  reduced  wave 

equation  which  is  obtained  from  (2.28)  by  substituting  F  (r)/^r  in 

n 

place  of  R  (r).  The  reduced  wave  function  F  (r)  satisfies 
n  n 


FJr)  ♦ 


(k  2(r)  +  -t-U  (r)  =  -  ^  (A  F  +  2B  F  )  .  (2.31) 

\  m  2  f  m  ZmmJ  mn  n  mn  n  v  ' 


The  radial  wave  functions  of  Eqs.  (2.28)  and  (2.31) 
satisfy  a  set  of  coupled  differential  equations.  The  exchange  of 
energy  between  normal  modes  or  the  mode  coupling  process  is 
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t! 

described  by  the  coupled  differential  equations  which  the  radial 

functions  satisfy.  The  RHS  of  Eqs.  (2.28)  and  (2.31)  can  be  thought 

of  as  a  spatial  distribution  of  sources  feeding  or  absorbing 

amplitude  to  or  from  mode  m.  The  coupling  coefficients,  A  and 

mn 

B  ,  are  therefore  proportional  to  radial  rates  of  interchange  of 
amplitude  between  modes  m  and  n. 

An  intuitive  picture  of  sound  propagation  via 
coupled  normal  modes  is  as  follows.  In  the  range  stratified  case, 
the  proportion  of  energy  alloted  to  each  normal  mode  of  oscillation 
remains  constant  in  range.  In  other  words,  the  partitioning  of 
energy  among  the  modes  does  not  vary.  When  a  range  dependence  of 
the  medium  is  introduced,  an  exchange  of  energy  between  modes 
occurs  and  the  partitioning  of  energy  between  the  modes  changes 
with  range  from  the  source.  The  increase  or  decrease  of  amplitude 
of  the  modes  as  a  function  of  range  is  described  by  the  behavior  of 
the  amplitude  of  the  radial  function. 

At  this  point  it  is  interesting  to  note  that  the 
mathematical  formalism  of  coupled  mode  theory  reduces  to  range 
invariant  mode  theory  as  the  range  dependence  of  the  medium  goes  to 
zero.  As  the  medium  approaches  horizontally  stratified,  the  <t>n(z;r) 
and  k^Cr)  become  range  independent.  When  this  occurs,  the  coupling 
coefficients  defined  by  Eqs.  (2.29)  and  (2.30)  vanish  and  Eqs. 
(2.28)  and  (2.24)  reduce  to  the  same  ones  encountered  in  conven¬ 
tional  mode  theory. 
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2.  Per  ]' vat  ion  of  Expressions  for  the  Coupling  Coefficients 

The  coupling  coefficients  defined  in  Eqs,  (2.29)  and 

(2.30)  are  important  functions  in  the  coupled  mode  theory 

formalism.  The  radial  functions  cannot  be  determined  until  the 

coupling  coefficients  as  a  function  of  range  are  known.  The 

definition  of  the  coupling  coefficients  given  in  Eqs.  (2.29)  and 

(2.30)  are  not  convenient  for  numerical  computation  because  they 

entail  calculation  of  range  derivatives  of  the  t}>n(z;r)  over  the 

entire  range  of  z  occurring  in  the  problem.  In  the  following 

analysis  expressions  for  A  and  B  more  amenable  to  numerical 
r  mn  mn 

calculation  are  derived. 

To  derive  an  alternate  expression  for  B  ,  one 

mn’ 

starts  with  Eq .  (2.24).  After  differentiating  with  respect  to 
range,  (2.24)  becomes 

2 

^2  i  +  2(kk  -  W  K  +  (k2  -  kn2)  ^n  =  °  •  (2-32) 


If  (2.32)  is  multiplied  by  p<t>m(z;r)  and  integrated  over  depth,  the 
result  is 


/ 


P<J>  — s-  <{>  dz  + 
m  gz2  Tn 


/ 


(k2~k^)pc)><J  dz 

n  m  n 


2k  k  6 
n  n  n,m 


/CO  # 

2pkk(f)  A 
mn 


(2.33) 

dz 


■  ■miu) mi  . . 


where  the  orthonormality  property  (Eq.  2.27)  of  the  <J>  has  been 
exploited.  To  illustrate  the  computational  technique,  a  two¬ 
layered  waveguide  geometry,  as  pictured  in  Fig.  II. 2,  will  be 
assumed.  The  extension  to  waveguides  having  more  layers  is 
straightforward.  Due  to  the  discontinuous  nature  of  the  depth 
functions  at  layer  interfaces,  the  integrals  from  0  to  °°  in  (2.33) 
are  shorthand  notation  for 


(2.34) 


With  Eq.  (2.34)  in  mind,  the  first  integral  in  Eq.  (2.33)  may  be 
integrated  twice  by  parts  with  respect  to  z  to  yield 


/ 


p*. 


dz 


f 


dz 


•  a  • 

P<t,n  T~rdZ  +  Pl4’lm9T*ln 
dz 


H(r) 


+  p2*2m  4  *2n 


lTln  dz  ylm 


H(r) 


H(r) 


(2.35) 


-  P 


2T2n  dz  y2m 


H(r) 


In  (2.35),  the  subscript  preceding  the  mode  index  on  the  $n  denotes 


FIGURE  1.2 

WAVEGUIDE  GEOMETRY  FOR  DERIVATION  OF  ALTERNATE 
EXPRESSIONS  FOR  COUPLING  COEFFICIENTS 
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the  depth  function  in  the  ith  layer.  The  boundary  conditions 
satisfied  by  4>n  are: 

'  0 

plh»  '  C2*2»|z.HW  (2-36) 

_a_ ,  =  i-  *  | 

3z  'lm  9z  2m  |  2=H(r-) 


Equation  (2.36)  combines  with  (2.35)  to  give 


/ 


,2  . 


P4> 


m  ~  2  n 
O  oz 


dz 


r 


T  il  *  plh»[3l  ♦in  '  3z  T2n 


fi-i.  -i-i. 


3z 


z=H(r) 


-  h  '  p2*2„] 


(2.37) 


z=H(r) 


Equation  (2.37)  can  be  used  in  (2.33)  to  get 


_  00 

[km2(r)  ‘  kn2(r)]  f0  P4>^n  dz  =  2knknVm  ‘  2Jo  dz 


r  a  •  a  :  1 

-  pl*lm[3l  *ln  '  9*  02nJz= 


H(r) 

(2.38) 


H(r) 


By  taking  the  total  range  derivative  of  Eqs.  (2.40)  and  (2.41)  it  is 
possible  to  derive  the  following  expressions. 
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-  (s,2V)kL. 


(2.43) 


Equations  (2.42)  and  (2.43)  combined  with  (2.39)  give 


B  (r) 
mn 


=  [kn2(r)  -  km2(r)]  1  1 2  /  Pkk<tm<t>n  dz  ♦  H(r)almn  J(2.44) 


where  a,  is  given  by 
lmn 


0.  (r) 

lmn 


3z  3z  l  1  J  z=H(r) 


(2.45) 


Equation  (2.44)  is  much  more  convenient  than  (2.30)  for  numerical 
calculation  because  it  contains  no  terms  involving  radial  deriva¬ 
tives  of  normal  mode  functions.  All  of  the  quantities  necessary  to 
evaluate  (2.44)  are  known  once  the  depth  functions  for  a  given  range 
have  been  computed. 

Equation  (2.44)  is  not  valid  for  m  =  n.  To  derive  an 

expression  for  B  ,  one  takes  the  range  derivative  of  the  mode 
r  nn 

orthogonality  relationship  as  follows. 


.  v1  MW«| 
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d  (Y,,(r)  r00 


p2*2.*2„  dz 


=  0 


(2.46) 


/  '  ’(Vlmhn  *  0l4’lr.iln)Js  *  /,  ,(p2*2n*2„  dz  *  P2*2m*2n)dZ 


H(r) 


=  -'Vlr/ln 


ll(r) 


+  H^2<l>2m^2n 


ll(r) 


The  integrals  in  Eq.  (2.46)  may  be  expressed  in  terms  of  the  B 
through  the  use  of  (2.30). 


mn 


B  (r)  ♦  B  (r)  =  2B  S  (r)  =  -lip, <fr}J>, n  (l-P./P.) 
mn  nm  mn  1  lm  In  1  / 


z-H(r) 


•  (2.47) 


In  (2.47)  the  superscript  S  means  the  symmetric  part  of  B  .  The 

mn 

diagonal  elements  of  B  are  easily  obtained  from  (2.47)  and  are 

mn 

given  by 


Hpij*ln] 

Bnn(r)  =  2  (1-Pj/P2) 


z=H(r) 


(2.48) 


Equation  (2.47)  can  also  be  derived  from  the  closure  relationship 
for  the  mode  functions.  Closure  is  a  property  of  complete,  ortho¬ 
normal  sets  of  functions  and  is  expressed  as 


7]  P<Pn(z; r)  4>n( z';r)  =  S(z-z') 


(2.49) 


z. 
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An  expression  for  A  can  be  obtained  by 

mn 

differentiating  the  integral  defining  with  respect  to  range. 


4—  B  (r) 

dr  mn 


■Hi 


H(r) 


P^Ir/ln  dz  + 


d=l 


(2.50) 


L(r)  =  /  PV*n  dz  +/  P^n  dz  +  A[Pl*lm*ln  '  p2*2,An]H 


(r) 


With  Eq.  (2.29),  (2.50)  becomes 


A  (r)  =  B  (r)  -  f  p4>  $  dz  -  Hp,<{>.  [<L  -  1  *(2.51) 

mnv  ;  mn1"  J  J  KlTlm[rln  2nJ|-i(r) 


It  is  possible  to  derive  an  expression  for  4>n  from  the  closure 
relation  of  Eq.  (2.49).  If  (2.49)  is  differentiated  with  respect  to 
range  multiplied  by<J>m(z;r)  and  then  integrated  over  z,  4>n  can  be 
shown  to  be  given  by 


VZ;r)  =  EniVz;r) 


(2.52) 


with 


Wr)  = 


+  Hpi^l2,^ln(1~Pl/p2^ 


H(r) 


(2.53) 
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In  deriving  Eq.  (2.52)  use  was  made  of  the  fact  that  the  density  has 
a  step  function  discontinuity  at  z=H(r),  i.e., 


p(z,r)  =  Pj  <D)  [H(r )  -  z]  -  P2©[z  -  H(r)] 


(2.54) 


( 


p(z,r)  =  H{p16(H(r)-z)  -  P2  6(z-ll(r)]  ’  (2.55) 


In  (2.54)  (H)(x)  is  the  Heavvside  unit  step  function  defined  by 


<Q)(x)  = 


1  x>0 


0  x<0 


Using  Eq.  (2.52)  in  (2.51),  one  may  obtain 


A  (r)  =  B  (r)  -  V  E  .  E  -  Hp.cJ),  (}.  -i 
mn  mn  nil  mil  1  lm  in 


2n 


.  (2.56) 


H(r) 


In  the  preceding  derivations,  any  range  dependence 
of  the  density  other  than  of  the  form  in  Eq.  (2.54)  was  not 
considered.  As  with  the  case  of  continuous  variation  of  density 
with  depth,  the  continuous  radial  variation  of  density  within  a 
layer  is  expected  to  be  of  negligible  importance.  Should  one  desire 
to  include  radial  density  variations  within  a  layer,  the  wave 
equation  for  the  acoustic  field  is  changed  (see  Ref.  2.13)  and  the 
preceding  developments  and  derivations  must  be  altered. 


4i’ 
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C.  Implementation  of  Coupled  Mode  Theory 

The  application  of  coupled  mode  theory  to  underwater 
sound  propagation  is  accomplished  in  four  stages.  The  first  stage 
entails  computation  of  sets  of  normal  mode  functions  on  a  mesh  of 
range  values  covering  the  radial  distance  of  interest.  In  the 
second  stage,  the  normal  mode  functions  are  used  to  compute  the 
coupling  coefficients  as  a  function  of  range  through  the  use  of  Eqs. 
(2.44)  and  (2.56).  In  the  third  step,  the  radial  functions  defined 
by  Eqs.  (2.28)  or  (2.31)  are  computed.  The  last  stage  combines  the 
normal  mode  and  radial  functions  according  to  Eq.  (2.23)  to  produce 
the  field.  In  each  of  stages  1  through  4  there  are  questions 
concerning  choice  of  the  range  mesh  appropriate  for  a  given 
problem.  These  questions  are  essentially  numerical  analysis 
concerns  and  are  not  considered  in  this  dissertation.  Suffice  it  to 
say,  the  range  sampling  must  be  fine  enough  to  allow  sufficiently 
accurate  numerical  solutions  of  the  differential  equations  for  the 
radial  functions. 

1 .  Solution  of  Depth  Equations 

The  numerical  solution  of  the  depth  equation  (2.24) 

subject  to  the  appropriate  boundary  conditions  is  a  numerical 

problem  familiar  to  many  areas  of  engineering  and  physics.  In 
2.14  2.15 

quantum  mechanics,  '  ’  *  for  example,  the  one-dimensional 

Schroedinger  equation  has  the  same  form  as  Eq.  (2.24)  where 
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and 


,  2.  2m 

k  (r’z)  ■  -^2  V(SJ 


.  2.  .  2m  c 

k„  <r)  =  rr  En 

n 


In  the  above  two  expressions  V(x)  is  the  potential  energy  felt  by  a 
particle  of  mass  m,  En  is  the  energy  eigenvalue,  and  "R  is  Plank's 
constant  divided  by  2tt.  Returning  now  to  underwater  acoustics,  the 
numerical  solution  of  Eq.  (2.24)  subject  to  the  previously 
described  boundary  conditions  is  a  familiar  and  well  documented 
procedure  and  will  not  be  discussed  in  any  detail.  References  2.2, 
2.7,  and  2.16  describe  three  different  methods  of  solution  of  Eq. 
(2.24). 


2 .  Solution  of  Radial  Equations 

The  solution  of  the  radial  equations  (2.28)  or 
(2.31)  is  not  a  routine  matter  and  will  be  discussed  in  more  detail 
here.  Assuming  that  the  coupling  coefficient  matrices  as  a 
function  of  range  are  known,  one  can  proceed  toward  a  solution  of 
the  radial  equations  in  a  number  of  ways.  A  method  used  in  making 
the  numerical  calculations  reported  later  in  this  dissertation  is 
an  iterative  technique.  In  this  approach,  the  set  of  coupled 
differential  equations  is  replaced  by  the  following  inhomogeneous 
equations  which  are  solved  iteratively. 


k  (r)  ♦ 
m 


,r2) 


(o)  _ 


=  0 


(2.57) 


y 
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p(i+l) 

m 


'  2 

k Cr) 

i  m 


4r“  m 


+  2B  F(lM*  (2.58) 

DI  ID  / 


The  iteration  process  on  Eq.  (2.58)  is  continued  until  the  pro¬ 
cedure  converges,  i.e.,  until  F^+^=F  ^  to  the  desired  accuracy. 
For  typical  range  variations  found  in  the  ocean,  the  coupling 

coefficients  are  expected  to  be  relatively  small  and  slowly 
varying.  For  this  reason,  it  is  expected  that  the  iterative 
procedure  just  described  will  converge  fairly  rapidly.  However,  as 
with  all  iterative  computational  schemes  one  must  be  prepared  to 
encounter  unstable  systems  of  equations  which  will  not  converge 
under  the  iterative  procedure. 

The  iterative  solution  of  the  radial  equations  is 
very  similar  to  the  "self-consistent  fields"  method  employed  in 
atomic  physics  to  compute  the  Hartree-Fock  atomic  wave  functions 
(see  Refs.  2.17  and  2.18).  The  Hartree-Fock  wave  functions  satisfy 
a  set  of  coupled,  integro-differential  equations;  however,  the 
iterative  process  employed  in  their  solution  is  identical  to  the 
one  just  described. 

A  second  type  of  iterative  procedure  that  is  very 
similar  to  the  one  described  in  Eqs.  (2.57)  and  (2.58)  is  a  Green's 
function  method.  Using  this  approach,  one  replaces  Eq.  (2.31)  with 
the  following  equation  which  is  solved  iteratively 


F(l  +  1-*(r)  = 
m 


^  CO 

Zf 

n  o  ' 


A  F(l) 
mn  n 


2B  F 
mn  m 


l)) 


(r.r')dr'  • 

(2.59) 


7,7 


The  Green's  function  g^°^  is  obtained  from  the  homogeneous 
equation  for  (Eq.  2.57)  and  is  given  by  (see  Ref.  2.10,  Ch.  5) 


E»°1(r,r'1  *  < 


Fl.frl  W’’ 

A(r') 


F2m^ 


r  >  r' 


(2.60) 


A(r') 


r  <  r' 


In  (2.60),  Fjm  and  F2m  are  two  linearly  independent  solutions  to 
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Eq.  (2.57).  This  technique  has  been  used  by  Williams  *  to 
calculate  radial  wave  functions  for  the  isovelocity,  wedge  shaped 
waveguide. 

Another  method  that  has  been  used  to  solve  coupled 
differential  equations  is  a  matrix  method.  In  this  technique,  a 
range  dependent  similarity  transformation  is  obtained  such  that 
Eqs.  (2.28)  or  (2.31)  are  diagonalized.  The  ensuing  system  of 
differential  equations  are  not  coupled  and  may  be  solved  by  conven¬ 
tional  methods.  After  solution  of  the  noncoupled  system  of 

equations,  the  desired  solutions,  i.e.,  the  R  or  F  ,  are  obtained 

m  m 

by  applying  the  inverse  similarity  transformation  to  the  noncoupled 
solutions.  More  details  concerning  this  method  may  be  obtained 
from  Refs.  2.20,  2.21,  and  2.22. 
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3.  Initial  Conditions 

Before  any  type  of  numerical  solution  for  the  radial 

functions  can  be  implemented,  the  proper  initial  conditions  of  the 

functions  must  be  known.  The  type  of  problem  in  underwater 

acoustics  most  likely  to  be  attacked  with  coupled  mode  theory 

involves  computation  of  the  field  due  to  a  monochromatic  point 

source  in  a  range  dependent  medium.  In  this  situation  it  can  be 

assumed  that  in  the  neighborhood  of  the  point  source,  the  field  is 

essentially  the  same  as  the  field  predicted  by  range  invariant  mode 

theory.  This  assumption  is  expected  to  be  quite  good  unless  severe 

range  variations  of  the  medium  occur  very  near  the  source. 

Excluding  these  pathological  cases  of  near-source,  severe  range 

variations,  one  obtains  the  initial  conditions  for  the  radial 

functions  by  requiring  that  they  match  the  field  for  a  range 

invariant  medium  at  some  point  r=r  which  is  close  to  the  source 

s 

region. 

Figure  II. 3  shows  a  portion  of  what  might  be  consid¬ 
ered  a  typical  range  dependent  waveguide.  In  the  notation  of  Fig. 
II. 3,  the  initial  conditions  for  the  radial  functions  and  their 
derivatives  are  taken  to  be  (see  Eqs.  2.21  and  2.23) 

R  (r  )  =  iTT  H^(k  r  )  <j)  (z  )  p(z  ) 
m  s  o  m  s'  m  o  v  o’ 

R  (r  )  =  p(z  )iir  k  H^(k  r  )  <h  (z  ) 
nr  s'  o'  mo  ^  m  s'  Tnr  o' 

(2.61) 

«-i*  km  HS1)(k»rs>  P(V>  • 


REGION  2 


SOURCE 

(0,  i„) 

* 


WATER  LAYER 


REGION  1 


SEDIMENT 


SUBSTRATE 


FIGURE  1.3 

WAVEGUIDE  STRUCTURE  USED  TO  DISCUSS  INITIAL  CONDITIONS  OF  RADIAL  EQUATIONS 
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In  Eq.  (2.61)  the  <f>  (z  )  is  a  mode  function  for  the  source  region 
mo 

(region  1)  computed  using  the  local  bathymetry  and  geoacoustic 

parameters.  With  the  initial  conditions  given  by  Eq.  (2.61),  the 

radial  functions  in  region  2  can  be  computed  numerically  using  one 

of  the  previously  mentioned  techniques.  Implicit  in  Eq.  (2.61)  is 

2  19 

the  neglect  of  backscattered  energy.  Previous  work  ’  has 
indicated  that  the  backscattered  portion  of  the  field  is  negligibly 
small.  For  this  reason  and  the  fact  that  the  inclusion  of  the 
backscattered  field  makes  a  formidable  numerical  problem  even  more 
difficult,  the  incoming  portion  of  the  acoustic'  field  generated  by 
the  range  dependence  of  the  medium  is  usually  neglected. 

4.  Attenuation  Effects 

Up  until  now,  no  mention  has  been  made  of  attenua¬ 
tion  effects.  What  is  customarily  done  in  conventional  mode  theory 
is  to  include  attenuation  effects  using  first  order  perturbation 
theory.  The  procedure  is  as  follows.  First,  it  is  assumed  that  the 
attenuation  effects  of  the  medium  may  be  treated  by  introducing  a 
small  imaginary  component  of  the  wave  number',  i.e.,  k(z)  becomes 


k(z)  +  ia(z) 


(2.62) 


a  <<  k 

where  a(z)  is  the  attenuation  as  a  function  of  depth  in  units  of 
nepers  per  unit  length.  Equation  (2.62)  is  next  used  in  (2.12)  to 


get 
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dz 


^  <}>n(z)  +  |^k 2 ( z )  -  k^2  +  2ia(z)k(z)j  <J>n(z) 


(2.63) 


where  a  has  been  neglected  as  a  second  order  quantity.  Conven- 

2.10  2.15 

tional  first  order  perturbation  theory  ‘  ’  *  can  be  applied  to 

2 

Eq.  (2.63)  to  determine  the  change  in  kR  due  to  the  perturbing  term 

2 

2ia(z)k(z).  The  change  in  k^  is  imaginary  and  is  given  by 


/*°°  2 
i<$  =  2i  /  ot(z) k(z)  p(z)  <t>  (z)  dz 

n  J  n 


(2.64) 


Hence,  the  new  eigenvalue  of  the  problem  becomes 


k  +  i6 


n  n 


_  ,  16 

,  2  1.  n 
k  I  1  +  — 


n 


(2.65) 


From  (2.65)  the  correction  to  k  can  be  found  and  is  seen  to  be 

n 


i6  \  i6 

k"C 


(2.66) 


Attenuation  effects  are  built  into  the  field  by 
retaining  the  unperturbed  mode  functions  while  using  Eq.  (2.66)  in 
the  radial  equations.  The  field  in  a  range  invariant  medium 
including  attenuation  effects  becomes  from  Eq.  (2.21) 
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*  ■  E  "o‘)(knr  *  TiT^V1’  *n(Jo> 


(2.67) 


If  one  uses  the  asymptotic  form  of  the  Hankel  function  valid  for 
k^rsl,  (2.67)  becomes 


IN  I 


-  6  r/2k 
n  n 


(2.68) 


n  A 


The  extension  of  this  perturbation  approach  to 

coupled  mode  theory  is  straightforward.  Since  in  coupled  mode 

2 

theory  the  eigenvalues  are  range  dependent,  so  too  are  the 

perturbation  corrections  to  k  (r).  Therefore  k  (r)  is  replaced 

n  n 


i6n(r) 

kn(r)  +  2kn(r) 


(2.69) 


in  analogy  with  Eq.  (2.66).  The  <5  ( r)  is  given  by 


/°° 

a(z,r)  k(z,r)  p(z)  4>n“(z;r)  dz  .  (2.70) 


Equation  (2.70)  is  simply  (2.64)  with  range  dependence  of  the 
modes,  wave  number,  and  attenuation  allowed  for.  When  the  depth 
functions  and  eigenvalues  for  a  given  range  have  been  computed,  the 
attenuation  coefficients  6  (r)  can  easily  be  calculated  via  (2.70). 
When  all  of  the  normal  mode  quantities  and  the  attenuation 


X. 


coefficients  have  been  determined  on  the  desired  range  grid,  the 
attenuation  effects  can  be  fed  into  the  radial  equations  (2.28)  or 
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(2.31)  by  replacing  k  (r)  with  Eq.  (2.69). 

ra 

In  the  preceding  developments  the  source  of  the 
attenuation  is  not  mentioned.  Whether  it  is  caused  by  thermo- 
viscous  effects,  relaxation  mechanisms,  or  any  other  number  of 
mechanisms  is  not  important.  One  must  simply  know  the  attenuation 
as  a  function  of  depth  and  range  for  the  waveguide  irregardless  of 
what  caused  it.  In  low  frequency  underwater  sound  propagation,  the 

attenuation  process  is  dominated  by  volume  absorption  effects  in 

•  .  2  23  f  24 

the  sediment  layers.  Published  experimental  work  *  is 

readily  available  documenting  the  results  of  measurement  of  a 

for  oceanic  sediments. 

5 .  Three-Dimensional  Coupled  Mode  Theory 

When  cylindrical  coordinates  and  azimuthal  symmetry 
are  not  appropriate  for  the  problem  at  hand,  a  more  general,  three- 
dimensional  coupled  mode  theory  of  the  following  form  can  be 
employed.  In  three-dimensional  cartesian  coordinates,  Eq.  (2.22) 
becomes 

V2'Kx,y,z)  +  k2(x,y,z)  ^(x.y.z)  =  0  •  (2.71) 

Employing  the  same  methods  as  used  in  section  P  of  this  chapter,  tjj 
is  assumed  to  be  given  by 


!'  .  - 
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iKx.y.z)  =  ^  Rn(x»y)  <i>n(z;x»y) 

n 


(2.72) 


with  R  and  *  satisfying 
n  Yn 


2 

3~R 

m 


'(x,y) 


R 

m 


+  B' 
mn 


V  R  \  (2.73) 
t  n  / 


r  ~2  7  2 

Aj-  ♦  k^(x,y,z)  -  kn  (x,y) 

,3z 


( z ;  x ,  y)  =  0  . 


(2.74) 


In  Eq.  (2.73)  V  is  the  transverse  gradient  operator  defined  by 


^3  ~  3 

Vt  =  X  37  +  y  ^  ’ 


(2.75) 


and  A'  and  B'  are  given  by 
mn  mn 


Amn(x’y)  =  f  P<t>mVt2<^n  d2 


(2.76) 


^  =  /  p4>mVnd2 


(2.77) 
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As  in  the  previous  case  the  4>^(z  ;x,y)  of  Eq.  (2.74)  satisfy 
continuity  of  p4>n  and  D4>n/9z  boundary  conditions  and  form  a  complete, 
orthonormal  set  of  functions. 

The  task  of  applying  coupled  mode  theory  to  a  truly 
three-dimensional  problem  is  at  best  a  formidable  one.  It  would 
entail  computation  of  .depth  functions,  eigenvalues,  and  coupling 
coefficient  matrices  over  a  two-dimensional  grid  of  (x,y)  points. 
The  amount  of  effort  that  would  be  required  to  accomplish  this  is 
probably  greater  than  would  be  expended  in  an  outright  three- 
dimensional  solution  of  Eq.  (2.71). 

One  type  of  problem  that  requires  a  cartesian 

coordinate  system  formulation  of  coupled  mode  theory  involves  a 

range  dependent  medium  in  which  the  lateral  variations  occur  in  one 

dimension  only.  This  type  of  waveguide  scenario  might  reasonably  be 

applied  to  acoustic  propagation  over  a  continental  slope  region  in 

which  the  bottom  slope  is  a  function  of  x  alone  (see  Fig.  II. 4).  To 

apply  coupled  mode  theory  to  problems  of  this  sort,  the  y- 

dependence ■ of  the  radial  equation  is  removed  by  taking  the  Fourier 

transform  of  Eq.  (2.73)  (see  Ref.  2.25).  The  Fourier  transform  of 

R  (x,y)  is  given  by 
in 

+oo 

Rm(x,y)  =J-J  dy  elYy  Rjx.y)  .  (2>78) 

m  v'Tn  ■'-«> 


Taking  the  transform  of  (2.73)  one  obtains 


SEA  SURFACE 


FIGURE  X4 

waveguide  geometry  for  medium  in  which  range  DEPENDENCE  OCCURS  IN  X- direction 

THIS  DEPICTION  IS  INVARIANT  WITH  RESPECT  TO  TRANSLATION  IN  THE  /-DIRECTION 
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(2.79) 


dx 


7  R  ( x ,  Y )  + 

i  m 


h>  -  -  £  (a;„  *  *;„ 


vB„ 


Equation  (2.79)  is  to  be  solved  over  a  sufficient  range  of  y  such 
that  the  inverse  transform  of  Rm(x,y)  can  be  computed  using 


R  (x,y)  =  —  f  dye'1Yy  R  (x,y)  .  (2.80) 

m  m  Joo  m 

Having  determined  Rm(x,y)  according  to  Eqs.  (2.79)  and  (2.80),  the 
field  is  constructed  from  Eq.  (2.72)  where  <t>n  now  depends  only  on  z 
and  x. 

At  this  point  it  is  worth  pointing  out  the 
differences  between  the  two-dimensional  coupled  mode  theory 
described  in  section  B  of  this  chapter  and  the  formalism  of  the 
theory  for  range  variations  occurring  in'  one  dimension  that  was 
just  described.  If  the  formalism  of  section  B  were  to  be  applied  to 
a  waveguide  as  pictured  in  Fig.  II. 4,  one  would  actually  be 
considering  a  bowl  shaped  structure  obtained  by  rotating  Fig.  II. 4 
about  an  axis  perpendicular  to  the  ocean  surface  and  passing 
through  the  source  position.  In  many  applications,  this  would  be  a 
poor  approximation  to  the  physical  system  being  modeled.  However, 
in  many  types  of  problems  in  which  the  source  is  some  distance 
removed  from  the  areas  of  significant  range  variation,  the 
differences  between  the  two  different  formalisms  of  coupled  mode 
theory  might  reasonably  be  expected  to  become  small.  The  reason 
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that  the  differences  are  expected  to  become  small  is  related  to  the 
fact  that  at  large  distances  from  the  source,  the  rotation  of  a 
waveguide  such  as  shown  in  Fig.  II. 4  about  the  source  position  is  a 
good  approximation  to  a  translation  in  the  y  direction,  even  for 
angles  of  rotation  that  are  fairly  large.  In  other  words,  the  bowl 
becomes  so  large  and  the  curvature  so  small  that  the  sound  field 
loses  knowledge  of  the  basin-like  structure  it  is  propagating  in 
and  becomes  a  good  approximation  to  the  field  in  the  waveguide  being 
modeled. 

The  advantage  of  applying  an  azimuthally  symmetric 
cylindrical  coordinate  system  whenever  possible  is  that  it  elimi¬ 
nates  the  need  to  perform  the  inverse  Fourier  transform  procedure 
of  Eq.  (2.80).  Since  the  inverse  transform  procedure  involves 
solution  of  (2.79)  for  many  values  of  y,  the  savings  in  computational 
effort  can  be  quite  significant.  A  distant  deep  water  source 
radiating  energy  toward  a  continental  slope  is  an  example  of  a 
propagation  problem  which  can  reasonably  be  handled  using  an 
azimuthally  symmetric,  cylindrical  coordinate  system  formulation. 


CHAPTER  III 

COUPLED  MODE  THEORY  IN  THE  ADIABATIC  APPROXIMATION 

In  this  chapter  the  adiabatic  approximation  to  coupled 
mode  theory  is  examined.  The  adiabatic  approximation  is  an 
approximation  on  the  set  of  coupled  radial  equations  in  which  all 
mode-mode  coupling  effects  are  ignored.  In  the  adiabatic  approxi¬ 
mation  the  amount  of  energy  initially  associated  with  a  mode 
remains  with  that  mode  since  energy  exchange  between  modes  has  been 
neglec  ted. 

The  adiabatic  approximation  to  coupled  mode  theory  has 

been  a  basic  part  of  the  formalism  of  the  theory  from  the  very 

31  .  3.2 

beginning.  In  the  two  benchmark  papers  by  Pierce  *  and  Milder 
in  which  the  full  framework  of  coupled  mode  theory  was  described,  it 
was  the  adiabatic  approximation  to  coupled  mode  theory  that  was 
actually  the  focal  point.  Pierce  applied  ray  theory  methods  to  the 
solution  of  the  adiabatic  range  equations  while  Milder  derived  ray 
theory  and  wave  theory  criteria  for  the  validity  of  the 
approximation.  A  large  portion  of  all  applications  of  coupled  mode 
theory  in  underwater  acoustics  has  been  within  the  adiabatic 
approximation  (see  Refs.  3.3-3. 7). 

The  material  contained  in  this  chapter  will  be  presented 
according  to  the  following  format.  First,  the  basics  of  the 
adiabatic  approximation  will  be  examined.  This  will  include  a 
discussion  of  energy  conservation  and  a  validity  criterion  for  the 
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approximation.  Second,  the  sensitivity  of  the  adiabatic  approxi¬ 
mation  to  the  range  dependence  of  the  ocean  bottom  will  be  examined. 
This  examination  will  take  place  within  the  context  of  a  realistic 
model  of  a  deep  water  oceanic  waveguide.  Finally,  the  field  in  the 
adiabatic  approximation  is  examined  with  regard  to  the  multipath 
conversion  processes  present  in  waveguides  having  nonhorizontal 
boundaries  and  interfaces.  This  will  involve  drawing  an  analogy 
between  wave  theory  and  ray  theory. 

A.  Background 

The  term  "adiabatic  approximation"  or  "adiabatic 

variation"  arises  in  many  branches  of  physics.  These  terms  were 

3  8-3.11 

originally  coined  in  the  Hamilton-Jacobi  '  ’  theory  of  action- 

angle  variables  in  classical  mechanics.  In  quantum  mechanics,  an 

adiabatic  approximation  (see  Ref.  3.12,  Chapter  8)  is  used  in  the 

solution  of  the  Schrodinger  equation  for  a  weakly  time-dependent 

Hamiltonian.  Another  quantum  mechanical  approximation  that  is 

mathematically  identical  to  the  adiabatic  approximation  of  coupled 
.  .  3  13-3  15 

mode  theory  is  the  Born-Oppenheimer  *  *  approximation  used  to 

separate  the  electron  and  nuclear  motions  in  the  solution  of 
Schrodinger 's  equation  in  atomic  physics  applications. 

In  the  application  to  coupled  mode  theory,  "adiabatic 
variation"  refers  to  the  range  dependence  of  the  properties  of  the 
acoustic  medium  in  question.  An  adiabatic  variation  of  some 
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property  of  a  waveguide,  e.g. ,  the  sound  speed,  is  a  finite  change 

in  that  parameter  occurring  such  that  locally 

8c  „ 

8r 

while  r  goes  to  infinity  in  a  manner  that 


=  Ac 


with  Ac  being  a  finite  change  in  the  parameter.  An  adiabatic 
variation  of  a  parameter  can  involve  quite  large  changes,  just  so 
long  as  the  local  rate  of  change  is  small.  The  "adiabatic 
approximation"  considered  in  this  chapter  is  based  upon  the 
assumption  of  adiabatic  range  variability  of  the  acoustic  medium 
and  is  expected  to  be  valid  in  the  limit  of  very  slow  local  radial 
rates  of  change  of  the  geoacoustic  parameters  of  the  waveguide. 


1 .  Fundamentals  of  the  Adiabatic  Approximation 

The  adiabatic  approximation  is  an  approximation 
applied  to  the  set  of  coupled  radial  equations  arising  from  the 
coupled  mode  theory.  The  approximation  is  based  on  the  assumption 
that  for  sufficiently  slow  range  variation  of  the  medium,  the 


effects  of  mode-mode  coupling  become  small  and  may  be  neglected. 
Based  on  this  assumption  the  radial  equations  in  the  adiabatic 
approximation  become  from  (2.28)  and  (2.31) 
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R  (r) 

R  (r)  ♦  -2— 


+  km  W  RmCr)  =  0 


(3.1) 


Fm(r)  +  (km2fr)  +  1/4r2)  Fm(r)  =  0 


(3.2) 


In  Eqs.  (3.1)  and  (3.2)  the  coupling  coefficients 

which  describe  the  mode-mode  interactions  have  been  neglected. 

From  Eqs.  (2.44)  and  (2.56)  it  is  apparent  that  the  B  and  A  are 

mn  mn 

proportional  to  first  and  second  range  derivatives  of  the 
properties  of  the  propagation  medium  (k(z,r)  and  H(r)).  Hence,  the 
neglect  of  the  coupling  coefficients  is  based  on  the  assumption 
that,  locally,  the  properties  of  the  medium  are  slowly  varying  and 
this  is  consistent  with  the  definition  of  adiabatic  variation. 

The  acoustic  field  computed  in  the  adiabatic 
approximation  differs  from  the  field  computed  using  the  full 
coupled  mode  formalism  only  in  the  radial  equations  ((3.1)  and 
(3.2)  versus  (2.28)  and  (2.31))  that  are  employed.  Of  course  the 
main  advantage  of  the  adiabatic  approximation  is  that  a  set  of 
uncoupled  radial  equations  is  to  be  solved  rather  than  a  coupled 
set.  Since  solving  uncoupled  equations  is  a  much  simpler  task  than 
solving  coupled  equations,  it  is  clearly  important  to  be  able  to 
take  advantage  of  the  adiabatic  approximation  whenever  possible.  A 
comparison  of  Eqs.  (2.5)  and  (3.1)  shows  that  the  adiabatic  range 
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equations  differ  from  the  horizontally  stratified  equations  only 
through  the  range  dependence  of  the  eigenvalue  k^Cr).  Therefore, 
it  is  expected  that  for  waveguides  having  slow  range  dependence, 
Eq.  (3.1)  will  exhibit  behavior  similar  to  a  Hankel  function. 

In  the  developments  of  Chapter  II  it  was  revealed 
that  coupled  mode  theory  when  applied  to  waveguides  having  range 
variable  boundaries  involves  an  approximation  to  the  physically 
proper  boundary  conditions  satisfied  by  the  field.  In  Chapter  IV  it 
will  be  shown  that  this  approximation  is  not  consistent  with  the 
theory  obtained.  In  the  adiabatic  approximation  to  coupled  mode 
theory  no  such  consistency  problem  arises  because  the  assumption  of 
adiabatic  variation  of  the  medium  is  of  the  same  order,  or  is 
consistent  with  the  approximation  to  the  physically  proper  boundary 
c  ond i t i on  s . 

The  partial  separation  of  radial  and  depth  variables 
effected  in  coupled  mode  theory  requires  that  the  normal  derivative 
of  the  field  at  a  nonhorizontal  boundary  be  approximated  with  the  z- 
derivative.  This  amounts  to  a  zero  order  approximation  of  Eq. 
(2.25). 


V  -  0(H)  .  (3.3) 

dZ 

.  *  ‘  2 

Since  the  coupling  coefficients  which  are  of  order  H  and  H  are 
neglected  in  the  adiabatic  approximation,  the  boundary  condition 
approximation  of  Eq.  (3.3)  is  consistent  with  the  assumption  of 


adiabatic  variation  of  the  medium.  The  significance  of  the 
consistency  will  become  more  apparent  in  Chapter  IV. 

In  the  first  part  of  this  chapter  it  was  stated  that 
the  energy  initially  associated  with  a  mode  remains  with  that  mode 
in  the  adiabatic  approximation.  To  illustrate  this  point  consider 
the  radial  energy  flux  f^  given  by  (see  Refs.  3.16  and  3.17) 

fr  =  i  Re  (pv*)  .  (3.4) 

In  (3.4)  p  is  the  acoustic  pressure  and  v  is  the  radial  component 
of  particle  velocity  given  by 

P  =  -  pff-  =  ia)P^  •  (3.5) 


vr  =  3?  •  (3.6) 

-icot 

In  (3.5)  the  time  dependence  of  the  field  was  taken  to  be  e  and 
"Re"  in  (3.4)  means  real  part  of.  For  the  field  expressed  in  the 
form 

'Kr , z)  =  5  Rn(r)  *nU:r)  ’  (3-7) 


Eq.  (3.4)  becomes 


where  the  dot  symbol  denotes  differentiation  with  respect  to  range. 
The  time  rate  of  change  of  energy  in  the  radial  direction,  Jr,  is 
obtained  by  integrating  (3.8)  over  the  surface  of  a  cylinder 
centered  about  the  depth  axis  at  range  r=0.  (See  Fig.  III.l.)  This 
yields 


j  =  f  2ttt  f 
r  /  r 


>r  ■  (\k  -  *v.)  *  g  (*«**.  ■  R*nR») '"»}  «• 


w>iere  the  definition  for  in  Eq.  (2.30)  along  with  (2.27)  were 


employed. 


In  the  adiabatic  approximation  all  quantities 


involving  coupling  coefficients  are  neglected;  hence,  (3.9)  becomes 
in  the  adiabatic  approximation 


:iwr/RJ*.  r*r  ) 
r  2  v  n  n  n  n  / 


(3.10) 


If  Eq.  (3.1)  is  multiplied  by  R*  and  subtracted  from  its  complex 
conjugate,  it  is  easy  to  show  that 


SOURCE 


c,U,  r),  p, 


CjU/  P2 


FIGURE  ELI 

WAVEGUIDE  GEOMETRY  FOR  DERIVATION  OF  RADIAL  ENERGY  FLUX 
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R  R* 
n  n 


R*R  + 
n  n 


•  • 

R  R*  -  R*R 
n  n  n  n 

r 


0 


d_ 

dr 


R  R* 
l  n  n 


1*R  ) 
n  n  / 


R  R*  -  R*R 
n  n  n  n 


r 


0 


(3.11) 


i[r(RR‘-R*R  1=0 
dr  [  \  n  n  n  n  /J 

Equation  (3.  II)  shows  that  in  the  adiabatic  approximation,  the 
quantity 


r 


'  R  R* 
k  n  n 


R*R 
n  n 


(3.12) 


is  constant  in  range  and  that  the  radial  flux  of  energy  through  a 
cylindrical  surface  transported  by  the  nth  mode,  \  is  range 
independent 


.(n)  i Trto  e 
r  =  ~2~  6n 


(3.13) 


This  verifies  that,  in  the  absence  of  attenuation,  the  amount  of 
energy  initially  associated  with  a  particular  mode  remains  with 
that  mode  in  the  adiabatic  approximation  to  coupled  mode  theory. 
Moreover,  the  total  radial  energy  flux  through  a  cylindrical 
surface  given  by 


J 


r 


(3.14) 


is  also  seen  to  be  range  independent,  indicating  that  energy  is 
conserved  within  the  context  of  adiabatic  mode  theory. 


I 


An  estimate  of  the  adiabatic  power  amplitude 
associated  with  a  given  normal  mode  may  be  made  based  on  the  WKB 
approximate  solution  of  Eq.  (3.1)  (see  Refs.  3.1,  3.12,  3.17,  3.18) 
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given  by 


i  f  k_(r')dr'  -  iir/4 


Rn(r)  -  i/2TT  p(zQ)  4>n(zo;0^> 


(3.15) 


/k  (r)r 
n 


In  Eq.  (3.15)  the  multiplicative  constants  were  chosen  such  that 
the  field  reduces  to  the  field  for  a  range  invariant  medium  near  the 
source  located  at  (zq,0)  (see  Chapter  II,  section  C-3).  Equations 
(3.12)— (3. 15 )  give 

=  2tt  0)P(z  )  <t>  (z  ;0)  •  (3.16) 

x*  OIlO 

2  2  2 

Jr  »  2tt  up(zo)  2  (Zo;0)  •  (3>17) 


Equations  (3.16)  and  (3.17)  are  the  adiabatic  power  amplitudes 
evaluated  within  the  WKB  approximation  and  depend  on  the  amplitude 
of  the  mode  functions  at  the  source  location. 
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2 •  Derivation  of  Adiabacity  Criterion 

In  this  subsection  of  Chapter  III  a  validity 
criterion  for  the  adiabatic  approximation  that  will  be  used  in 
later  sections  of  this  chapter  is  derived.  This  criterion  is  based 

on  coupled  mode  wave  theory  and  was  originally  derived  by 

3  2  . 

Milder.' ’  The  derivation  presented  here  is  a  reproduction  of 

Milder's  derivation.  Other  validity  criteria  for  adiabacity  have 

been  derived  (see  Refs.  3.4  and  3.6);  however,  Milder's  criterion 

is  of  a  form  most  suitable  for  the  investigations  presented  later  in 

this  chapter. 

The  derivation  starts  with  Eq.  (2.31) 

approximated  as 

+  k  2 (r)l  F  (r)  =  -  2  V  B  (r)  F  (r)  •  (3.18) 

,2  m  v  1  \  m  mn  n 

L  dr  ->  n 


In  writing  (2.31)  as  (3.18)  it  has  been  assumed  that  the  range 

2 

variations  of  interest  are  far  enough  from  the  source  that  the  l/4r 
term  may  be  neglected  and  that  Affln  which  is  second  order  in  the 
range  dependence  of  the  medium  is  negligible  compared  to  B  .  The 
solution  to  Eq.  (3.18)  can  be  expressed  in  terms  of  a  homogeneous 
solution,  and  a  particular  solution, 


Fm^ 


F(0)  +  Ffl) 

m  m 


(3.19) 
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The  solution  is  of  course  the  adiabatic  solution.  The 

m 

particular  solution  can  be  written  in  terms  of  a  Green's 

m 

function  as 


:CD 


Fi"<r>  ■  2Hf  fn(r'>*™(r'r'>dr' 

n  •'o 


(3.20) 


In  (3.20)  is  the  Green's  function  appropriate  to  the  homogeneous 
form  of  (3.18).  In  evaluating  (3.20),  WKB  approximations  for  the 
quantities  in  the  integral  are  employed.  These  are 


gm(r’r,)  =  2^TeXP[ijf,  km(x)dX] 


(3.21) 


Fn(r') 

-  exp|i  j 

-1"  1 

kn(x)dx| 

(3.22) 

1 

| 

Fn(r') 

~  ik  (r') 
n 

Fn(r')  . 

(3.23) 

Equations  (3.21)  -  (3.23) 

in  (3.20) 

give 

n  */o 

B  (r')ik 
mn  n 

""’e  \f 

[kn(x)  -  km(x)]dx} 

2ikm<r’) 

p  (  7 

x  exp  |  i  f  k 

m(x*dx}dr 

t 

(3.24) 

- 
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In  order  to  derive  a  criterion  that  does  not  depend 

on  the  particular  details  of  a  given  propagation  medium,  the 

following  assumptions  are  made.  First,  it  is  assumed  that  the  range 

variations  of  the  medium  occur  only  over  a  region  of  space  between  0 

and  r  .  Second,  it  is  assumed  that  the  eigenvalues  and  B  are 
o  mn 

slowly  varying  functions  of  range.  In  problems  involving  realistic 
deep  water,  oceanic  waveguides  this  assumption  is  quite  good. 
Under  the  above  two  assumptions  Eq.  (3.24)  becomes 


Fi1) (T) 


ik  r  r 

E  •  '  / 


o  i(k  -k  )r' 
e  dr' 


(3.25) 


..  i(k  -k  )r 

Lik  r  m  n  o  , 
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B  (r)  e 
mn 


ifk-’O 

m  n 


(3.26) 


In  obtaining  (3.26)  the  ratio  k  /k  in  (3.24),  which  is  a  term  of 

n  tn 

order  unity,  was  replaced  by  1. 

Milder's  criterion  is  obtained  from  (3.26)  by 
requiring  that  the  magnitude  of  be  much  less  than  the  magnitude 

of  the  WKB  form  of  adiabatic  solution  in  Eq.  (3.22).  Hence,  for  the 
effects  of  mode  coupling  to  be  small, 


Fm1)(r) 


«  1 


(3.27) 


Equation  (3.27)  will  be  satisfied  if,  term  by  term  from  (3.26) 
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2B  (r) 
mn 


k  (r)  -  k  (r) 
m  n 


«  1 


(3.28) 


Equation  (3.28)  states  that  if  the  coupling  coefficient  Bmn 
describing  the  mode  interaction  between  modes  m  and  n  is  small 
compared  to  the  difference  in  eigenvalues,  then  adiabatic  mode 
invariance  is  assured.  Equation  (3.28)  is  to  be  satisfied  for  all 
combinations  of  m  and  n,  especially  for  those  where  m  and  n  differ 
by  1  and  k  -k  is  the  least.  A  somewhat  relaxed  version  of  (3.28) 
is  obtained  by  requiring  the  coupling  between  nearest  neighbor 
modes  to  be  small. 


2B  . 
m,m+l 

k  -  k  . 
m  m+1 


«  1 


(3.29) 


Equation  ( 3 . 29")  is  the  adiabacity  criterion  that 
will  be  used  in  the  next  section  of  this  chapter.  In  terms  of  the 
mode  cycle  distance,  xm,  defined  by 


x 


m 


2-rr 

k  -  k  . 
m  m+1 


(3.30) 


Eq.  (3.29)  becomes 


B  ,  x 
m,m+l  m 


=  £  «  1 


(3.31) 


where  e  is  a  smallness  parameter. 


63 


3.  Types  of  Problems  Suitable  for  the  Adiabatic 

Approximation 

The  requirement  of  adiabatic  range  dependence  of  the 
medium  narrows  the  scope  of  problems  that  may  be  attacked  using 
adiabatic  mode  theory.  Since  interchange  of  energy  is  not  allowed 
in  adiabatic  mode  theory,  it  is  probably  not  suited  for  use  in 
waveguides  in  which  the  number  of  propagating  modes  varies 
drastically  in  range.  Such  a  situation  would  occur  if  adiabatic 
mode  theory  were  applied  to  a  continental  slope  propagation 
problem. 

For  a  20  Hz  source  located  in  deep  water  (5000  m) 
there  would  typically  be  about  110  propagating  modes.  For  a 
receiver  located  on  the  continental  shelf  in  shallow  water  (700  m) 
there  would  exist  about  20  propagating  modes.  In  this  situation  one 
could  not  expect  the  adiabatic  approximation  to  work  very  well  over 
extended  ranges  since  the  energy  associated  with  a  particular  mode 
disappears  from  the  problem  when  that  mode  reaches  cutoff.  Clearly 
the  process  of  mode  conversion  in  which  energy  is  transferred  among 
the  modes  would  be  important  in  this  case. 

The  types  of  problems  that  are  likely  candidates  for 

the  use  of  adiabatic  mode  theory  are  those  in  which  the  principle  of 

conservation  of  mode  index  applies.  The  principle  of  conservation 

3.1 

of  mode  index  was  used  by  Pierce  ’  to  describe  propagation 
situations  in  which  the  number  of  propagating  modes  remained 
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constant  in  range.  In  this  case,  the  range  dependence  of  the  medium 
is  assumed  to  be  accommodated  by  small  adjustments  of  the  mode 
functions  in  range  without  the  need  of  energy  transfer. 

Conservation  of  mode  index  can  reasonably  be  expected  to  apply  to 
waveguides  having  relatively  horizontal  boundaries  and  geoacoustic 
parameters  that  are  weakly  dependent  on  range.  Since  the  number  of 
propagating  modes  is  particularly  sensitive  to  water  depth,  wave¬ 
guides  having  large  depth  variations  cannot  be  treated  as 

adiabatic.  On  the  other  hand,  if  one  wanted  to  describe  the 

propagation  of  sound  over  a  relatively  flat  ocean  basin  where  the 
sediment  composition  varied  slowly  with  range,  then  the  adiabatic 
approximation  is  likely  to  be  quite  good. 

B.  Sensitivity  of  Adiabatic  Approximation  to  Ocean  Bottom 

Properties 

In  the  study  of  underwater  sound  propagation  in  the  ocean 
it  is  customary  to  consider  the  waveguide  in  terms  of  two  component 
parts,  the  water  layer  and  the  bottom.  The  bottom  is  considered  to 
begin  at  the  water-sediment  interface  and  continue  into  the 

substrate  to  infinity.  Depending  on  the  frequency  of  interest,  the 
bottom  can  influence  the  field  in  the  water  column  to  a  great 
extent.  At  low  frequencies  where  normal  mode  techniques  are  most 
applicable,  the  bottom  is  an  important  part  of  the  propagation 


medium. 
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In  this  section  of  Chapter  III  the  sensitivity  of  the 
adiabatic  approximation  to  the  range  dependence  of  the  ocean  bottom 
is  examined.  Previous  work  along  these  lines  (see  Refs.  3.2  -  3.5) 
has  focused  on  range  variability  in  the  water  column  alone.  In  most 
realistic,  deep  water  oceanic  waveguides  the  water  layer  is 
sufficiently  slowly  varying  with  range  that  the  adiabatic  approxi¬ 
mation  is  valid.  Milder  showed  that  for  typical  horizontal  sound 
speed  gradients  in  the  ocean,  adiabatic  propagation  in  the  water 
layer  was  assured  for  frequencies  less  than  600  Hz. 

The  ocean  bottom  is  not  nearly  as  stable  as  the  water 
column  and  much  larger  range  variability  is  present.  The  work 
presented  here  will  focus  on  two  range  dependent  bottom  parameters, 
sediment  sound  speed  and  water-sediment  interface  slope,  and  will 
explore  how  the  adiabatic  approximation  depends  on  these 
parameters. 

The  method  used  to  examine  the  sensitivity  of  the 

adiabatic  approximation  to  the  bottom  parameters  is  through  the 

use  of  Milder's  criterion  given  in  Eq.  (3.31).  The  procedure  is  to 

select  a  value  for  £  in  (3.31)  and  compute,  for  each  mode  number  m, 

the  magnitude  of  B  .  at  which  the  equality  in  (3.31)  is 

m*  tn+ 1 

satisfied.  The  magnitude  of  B  .  depends  on  the  properties  of  the 

m,m+i 

bottom  in  question.  This  allows  one  to  determine  the  values  of 
sediment  radial  sound  speed  gradient  and  water-sediment  interface 
slope  that  may  be  tolerated  within  the  adiabatic  approximation. 


In  the  choice  of  e  in  (3.31)  one  is  allowed  a  certain 
amount  of  latitude.  Only  when  the  adiabatic  solution  can  be 
directly  compared  with  the  full  coupled  mode  solution  will  it  be 
possible  to  determine  the  proper  value  of  the  smallness 
parameter  £  to  ensure  adiabacity.  Since  a  coupled  mode  theory 
calculation  for  a  general  oceanic  waveguide  is  not  yet  able  to  be 
done,  it  is  possible  only  to  select  what  is  considered  to  be  a  safe 
value  for  E  and  proceed  with  the  investigation.  Hence,  it  is  not 
now  possible  to  determine  absolutely  the  magnitude  of  bottom  slope 
and  radial  sound  speed  gradient  that  may  be  treated  adiabatically 
though  the  relative  sensitivity  of  the  adiabatic  approximation 
among  different  bottom  types  can  be  studied  and  understood.  The 
value  of  £  used  in  the  studies  presented  here  is  £=0.1. 


1 .  Dependence  on  Radial  Sound  Speed  Gradient 

The  first  range  dependent  parameter  to  be  investi¬ 
gated  is  the  sediment  sound  speed.  Figure  III. 2  depicts  the 
waveguide  geometry  used  in  the  numerical  calculations  that  were 
performed.  The  first  layer  is  water  4500  m  deep  having  an  arbitrary 
sound  speed  profile.  (For  the  calculations  presented  here,  a  mid- 
Pacific  type  profile  was  used. )  The  sediment  layer  is  300  m  deep 
and  is  assumed  to  depend  on  range  and  depth  as 
c2(r,z)  =  cQ(r)  +  g(z-H) 

In  other  words,  the  sound  speed  varies  with  range  and  depth  in  such 
a  way  that  the  depth  variation  is  linear  and  the  radial  gradient  of 
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sound  speed  is  depth  independent.  The  sediment  layer  will  be 

treated  as  clay,  silt,  and  sand.  These  sediment  types  span  the 

spectrum  of  sediment  porosities  usually  found  in  oceanic  sediments. 

The  substrate  is  treated  as  an  infinite  isovelocity  liquid  half- 

3.19 

space.  The  geoacoustic  parameters  for  the  layers  in  Fig.  III. 2 

are  listed  in  Table  I. 

Figure  III. 3  depicts  the  relationship  between  mode 
number  and  local  radial  sound  speed  derivative  such  that  Eq.  (3.31) 
is  satisfied  for  e=0.1.  The  frequency  is  20  Hz  and  the  sediment  is 
treated  as  a  clay.  There  are  106  discreet  modes  for  this  frequency 
and  layering  geometry.  Two  curves  are  present  in  Fig.  III. 3.  The 
solid  curve  corresponds  to  no  depth  dependence  of  the  sediment 
sound  speed,  i.e.,  isovelocity  in  depth,  while  the  dotted  curve 
corresponds  to  a  constant  depth  gradient  of  1  sec 

The  curves  in  Fig.  III. 3  were  generated  from  Eq. 
(2.44)  using  H=0. 

'  rk  * !  T7]  l  2oki;Wz 

Since  the  range  variations  are  restricted  to  the  sediment  layer  and 
are  given  by 


c2(r,z)  =  cQ(r)  +  g(z-H) 


(3.33) 


Eq.  (3.32)  becomes 


TABLE  I 


SEDIMENT  AND  SUBSTRATE  GEOACOUSTIC  PARAMETERS 


Bottom  Type 

Sound  Speed 
(m/sec) 

Density 

(g/cc) 

Clay 

1526.6 

1.27 

Silt 

1563.0 

1.60 

Sand 

1723.8 

2.08 

Substrate 

2700.0 

2.60 

CLAY,  20  Hz 


MAY  BE  TOLERATED  WITHIN  THE  ADIABATIC  APPROXIMATION 
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B(r)  = 
mn 


-2uj  cQ(r)p2 


/> 


(z;r)  <t>(z;r) 


-  k  “1  sed 

m  J 


dz 


(3.34) 


c2  (r,z) 


From  Eq.  (3.34)  it  is  straightforward  to  determine  the  values  of 

CQ(r)  that  cause  Eq.  (3.31)  to  be  satisfied  for  a  particular  choice 

of  £.  In  this  manner  Fig.  III. 3  and  all  similar  curves  were 

generated.  The  mode  functions,  (j)^,  used  in  the  computation  of 

3.20  3.21 

(3.34)  were  obtained  by  numerical  integration  ’  of  the  depth 
equation  (2.24). 

There  are  several  interesting  features  in  Fig. 
III. 3.  First,  note  the  large  dip  in  the  curve  for  g=0  at  about  mode 
number  24.  This  dip  marks  the  dividing  line  between  modes  with 
turning  points  in  the  water  and  those  with  turning  points  in  the 
sediment.  The  term  "turning  point"  refers  to  points  at  which  the 

2  2  2 

vertical  wavenumber  K  (z;r)=k  (r,z)-k  (r)  in  (3.35)  becomes  zero. 

n  n 


[^2  +  Kn(z;r)]*n(z; 


;r)  =  0 


(3.35) 


Hence,  a  turning  point  separates  regions  in  which  the  mode 
functions  4>n  are  oscillatory  and  exponential  in  behavior.  For  a 
typical  deep  water  sound  speed  profile  as  shown  in  Fig.  III. 2,  the 
modes  of  low  index  number  have  turning  points  wholly  within  the 
water  column  and  are  exponential  in  the  bottom  region.  These  modes 
are  said  to  be  trapped  in  the  water  column.  As  mode  index 
increases,  the  lower  turning  points  begin  to  occur  within  the 
bottom  and  oscillatory  behavior  of  the  mode  function  occurs  over  a 
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portion  of  the  bottom.  Such  modes  are  said  to  be  bottom  interacting 
modes. 

In  Fig.  III. 3,  the  modes  of  index  less  than  24  are 
trapped  in  the  water  columns  and  are  exponential  in  the  bottom 
region  where  the  range  changes  are  occurring.  Hence,  these  modes 
are  insensitive  to  what  is  happening  in  the  bottom  and  may  be 
treated  adiabatically  for  quite  large  radial  sound  speed  gradients. 

Another  feature  of  Fig.  III. 3  is  that  the  presence 
of  a  depth  gradient  washes  out  the  large  dip  beginning  at  mode 
number  24  and  approaches  the  curve  for  g=0  at  the  higher  mode 
numbers.  The  general  features  of  both  curves  in  Fig.  III. 3  may  be 
understood  with  reference  to  Fig.  III. 4.  When  a  mode  has  a  turning 
point  in  the  sediment,  it  may  be  viewed  as  a  ray  of  equivalent 
bottom  angle,  0^,  given  by 


0 

m 


-1 

cos 


CjCH) 


V 

m 


V  =  u)/k 
m  m 


where  is  the  phase  velocity  of  the  nth  mode.  The  equivalent 
bottom  grazing  angle  is  an  increasing  function  of  mode  number. 

The  large  dip  in  the  solid  curve  of  Fig.  III. 3  at 
mode  number  24  corresponds  to  the  first  mode  with  a  turning  point  in 
the  sediment  layer.  This  mode  corresponds  to  an  equivalent  ray  of 
very  shallow  grazing  angle.  When  no  depth  gradients  are  present  in 
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the  bottom,  this  ray  covers  the  largest  radial  distance  through  the 
sediment  on  its  path  to  the  substrate  and  back  to  the  water  column. 
This  equivalent  ray,  therefore,  is  exposed  to  the  largest  relative 
radial  change  in  sound  speed  as  it  cycles  through  the  sediment 
layer.  As  the  mode  number  increases,  the  equivalent  ray  bottom 
angle  increases  and  the  radial  distance  covered  by  the  ray  as  it 
cycles  through  the  sediment  layer  decreases  (see  top  row  in  Fig. 
III. 4).  The  adiabatic  approximation  depends  on  the  relative  amount 
of  radial  change  in  sediment  sound  speed  over  the  mode  cycle 
distance.  Hence,  for  the  case  of  no  vertical  sound  speed  gradient, 
the  higher  modes  should  be  able  to  withstand  larger  radial  sound 
speed  gradients  and  remain  adiabatic.  This  is  precisely  the 
behavior  of  the  solid  curve  in  Fig.  III. 3.  One  sees  that  the  amount 
of  radial  sound  speed  gradient  that  may  be  tolerated  in  the 
adiabatic  approximation  increases  steadily  with  mode  number. 

When  a  vertical  sound  speed  gradient  is  present  in 
the  bottom,  the  previous  situation  is  altered  somewhat  (see  bottom 
row  in  Fig.  III. 4).  When  g=0,  the  mode  with  the  shallowest 
equivalent  bottom  angle  no  longer  corresponds  to  the  largest  cycle 
distance  through  the  sediment  due  to  the  refraction  caused  by  the 
vertical  sound  speed  gradient.  As  the  mode  number  and  hence  the 
equivalent  bottom  angle  is  increased,  the  cycle  distance  through 
the  sediment  begins  to  increase  until  the  mode  turning  points  reach 
the  substrate.  At  this  time,  the  cycle  distances  begin  to  decrease 
and  approach  the  cycle  distances  for  g=0  as  the  equivalent  bottom 
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angles  become  so  steep  that  the  refractive  effects  of  the  vertical 
gradient  become  negligible.  This  is  the  behavior  exhibited  by  the 
dotted  curve  in  Fig.  III. 3.  The  large  dip  around  mode  24  vanishes 
due  to  the  decreased  cycle  distances  of  these  modes  when  a  vertical 
gradient  is  present  and  the  dotted  curve  approaches  the  solid  curve 
as  mode  number  increases  and  the  cycle  distances  for  the  two  cases 
become  similar. 

Figures  III. 5  and  III. 6  are  the  same  type  curves  as 
Fig.  III. 3  except  that  the  sediment  is  treated  as  silt  and  sand. 
Figures  III. 5  and  III. 6  exhibit  the  same  type  of  behavior  as  Fig. 
III. 3  except  that  the  large  dip  in  the  solid  curves  corresponding  to 
the  first  mode  with  a  turning  point  in  the  bottom  occurs  at  higher 
mode  numbers  because  of  higher  surficial  sound  speed  values  in  silt 
and  clay. 

Figure  III. 7  is  a  superposition  of  Figs.  III. 3, 
III. 5,  and  III. 6  for  g=l.  Based  on  Fig.  III. 7  and  a  choice 

of  £=0.1,  it  can  be  said  that  a  radial  sediment  sound  speed  gradient 

-3-1  . 

of  about  10  sec  could  be  tolerated  within  the  adiabatic 

approximation  for  all  the  sediment  types  and  all  mode  numbers. 

In  the  lower  mode  number  portions  (less  than  m  60) 
of  Fig.  III. 7,  the  curves  lie  in  the  order  sand,  silt,  and  clay.  In 
this  low  mode  region  the  sand  can  withstand  the  highest  radial 
gradients  followed  by  silt,  and  then  clay.  This  behavior  is  related 
to  the  fact  that  the  characteristic  acoustic  impedance  of  sand  is 


AND  MAXIMUM  RADIAL  SEDIMENT  SOUND  SPEED  DERIVATIVE  THAT 
MAY  BE  TOLERATED  WITHIN  THE  ADIABATIC  APPROXIMATION 


CURVE  REPRESENTING  THE  RELATIONSHIP  BETWEEN  MODE  NUMBER 
AND  MAXIMUM  RADIAL  SEDIMENT  SOUND  SPEED  DERIVATIVE  THAT 
MAY  BE  TOLERATED  WITHIN  THE  ADIABATIC  APPROXIMATION 
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higher  than  that  of  silt,  which  is  higher  than  that  of  clay.  Hence, 
for  the  low  mode  numbers  in  question,  the  sand  bottom  has  the  lowest 
portion  of  acoustic  energy  present  in  the  bottom  and  is  least 
affected  by  the  presence  of  a  radial  sediment  sound  speed  gradient. 
Silt  is  affected  to  a  greater  degree  since  its  characteristic 
impedance  is  less  than  that  of  sand.  Clay  has  the  smallest 
impedance  and  is  affected  to  the  greatest  degree. 

In  long  range  propagation  problems  it  is  likely  that 
only  the  lower  modes  are  of  importance  at  long  ranges  from  the 
source  because  of  attenuation  of  the  bottom  interacting  modes.  In  a 
situation  such  as  this,  the  magnitude  of  radial  sound  speed 
gradient  that  may  be  handled  within  the  adiabatic  approximation  can 
be  seen  to  vary  over  an  order  of  magnitude  as  the  bottom  type  varies 
from  clay  to  sand. 

2.  Dependence  on  Local  Bottom  Slope 

The  next  parameter  isolated  and  examined  is  the 
local  bottom  slope.  The  waveguide  geometry  that  was  used  for  the 
calculation  of  the  normal  modes,  coupling  coefficients,  and  mode 
cycle  distances  is  shown  in  Fig.  III. 8.  In  this  figure  there  is  a 
water  layer  of  variable  depth  given  by  H(r)  overlying  isovelocity 
sediment  and  substrate  layers.  The  calculations  to  be  presented 
were  performed  at  a  range  where  the  water  depth  was  4500  m  and  the 
sediment  depth  300  m.  Again,  the  water  sound  speed  profile  is  a 
typical  mid-Pacific  type  profile. 


SEDIMENT 

c2'  P 2 


dr 


SUBSTRATE 

c3'  P3 


FIGURE  111.8 

WAVEGUIDE  LAYER  AND  SOUND  SPEED  STRUCTURE 
FOR  EXAMINATION  OF  THE  EFFECTS  OF  A  SLOPING 
WATER-SEDIMENT  INTERFACE 
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Figure  III. 9  shows  the  relationship  between  mode 
number  and  local  bottom  slope,  dH/dr,  such  that  Eq.  (3.31)  is 
satisfied  for  £=0.1,  with  the  sediment  being  treated  as  a  clay.  For 
a  given  mode  number  m,  values  of  dH/dr  lying  beneath  the  curve 
correspond  to  slopes  that  may  be  treated  adiabatically  while  those 
lying  above  the  curve  correspond  to  slopes  that  may  not  be  treated 
adiabatically.  Based  on  the  choice  of  £=0.1,  Fig.  III. 9  reveals 
that  most  of  the  modes  propagate  adiabatically  for  slopes  up  to 
about  0.5  °. 


As  was  done  in  the  previous  case,  the  curve  in  Fig. 
III. 9  was  generated  from  Eq.  (2.44),  this  time  with  k=0. 

Bmn  '  ~  2  _  ”“2  H^°lmn  (3.36) 
n  m 

From  (3.36)  and  the  expression  for  o  given  in  (2.45),  the  values 
of  H  that  satisfy  the  equality  in  (3.31)  can  be  found  once  the  mode 
functions  and  eigenvalues  are  known. 

The  curves  in  Figs.  III. 10  and  III. 11  are  the  same 
type  as  the  one  shown  in  Fig.  III. 9  except  that  the  sediment  layer 
is  treated  as  silt  and  sand,  respectively.  Figures  III. 9,  III. 10, 
and  III. 11  are  similar  in  that  they  are  oscillatory  over  a  large 
portion  of  the  mode  numbers  present.  The  oscillatory  behavior 
begins  at  the  mode  number  corresponding  to  the  first  mode  with  a 
turning  point  in  the  sediment.  The  oscillatory  regions  of  the 


WITHIN  THE  ADIABATIC  APPROXIMATION 


CURVE  REPRESENTING  THE  RELATIONSHIP  BETWEEN  MODE  NUMBER 
AND  MAXIMUM  LOCAL  BOTTOM  SLOPE  THAT  MAY  BE  TOLERATED 
WITHIN  THE  ADIABATIC  APPROXIMATION 
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curves  for  silt  and  sand  are  delayed  in  mode  number  relative  to  the 
clay  case  because  of  critical  angle  effects  present  in  silt  and  sand 
that  increase  the  mode  number  at  which  bottom  penetration  first 
occurs. 

Figure  III. 12  is  a  superposition  of  Figs.  III. 9, 
III.  10,  and  III. 11.  Note  that  the  curves  tend  to  lie  in  the  order 
clay,  silt,  and  sand,  with  clay  being  able  to  tolerate  more  slope 
angle  than  sand.  This  ordering  follows  the  ordering  of  the 
characteristic  acoustic  impedances  of  the  sediment,  with  clay 
having  the  smallest  and  sand  the  largest.  This  behavior  is  not 
inconsistent  with  the  notion  that  an  acoustically  hard  sloping 
bottom  would  give  rise  to  more  mode  conversion  than  an  acoustically 
softer,  more  penetrable  bottom.  The  curves  in  Fig.  III. 12,  for  the 
choice  £=0.1,  indicate  that  mode  conversion  effects  begin  to  arise 

•  O  O 

in  sand  for  local  bottom  slopes  between  0.057  and  0.57  whereas 
for  clay,  mode  conversion  effects  begin  at  slope  angles  on  the  order 
of  0.5°  to  1.0°.  Hence,  the  mode  conversion  process  due  to  a 
sloping  bottom  seems  to  be  fairly  sensitive  to  bottom  type. 

3.  Summary 

The  dependence  of  the  adiabatic  approximation  on  the 
geoacoustic  parameters  of  the  bottom  has  been  examined  using  a 
criterion  developed  by  Milder  which  requires  that  the  amount  of 
mode  coupling  occurring  over  a  mode  cycle  distance  be  small.  It  has 
been  shown  that  the  magnitude  of  local  sound  speed  gradient  that  may 
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be  handled  within  the  adiabatic  approximation  can  vary  over  more 
than  an  order  of  magnitude  as  the  sediment  type  varies  from  clay  to 
sand.  It  has  also  been  shown  that  the  magnitude  of  local  bottom 
slope  can  vary  over  an  order  of  magnitude  as  the  sediment  changes 
from  clay  to  sand;  the  acoustically  softer  clay  bottom  causes  less 
mode  conversion  than  the  harder  sand  bottom. 

C.  Multipath  Information  Within  the  Adiabatic  Approximation 

3  92 

Recently,  Tindle  and  Guthrie  ’  published  a  paper 

exploring  the  relationship  between  ray  theory  and  normal  mode 

theory  in  horizontally  stratified  acoustic  propagation  geometries. 

Their  basic  conclusion  was  that  raypaths  correspond  to  the 

constructive  interference  of  a  group  of  normal  modes  centered  about 

the  mode  of  index  that  corresponds  to  the  equivalent  ray  angle  of 

the  ray  path  in  question.  This  conclusion  allowed  them  to  draw  a 

connection  between  ray  theory,  with  its  physically  intuitive  notion 

of  raypaths,  and  mode  theory,  which  is  exact  but  often  not  useful 

for  gaining  insight  into  the  physical  processes  which  influence  the 

propagation  of  underwater  sound.  The  idea  of  raypaths  arising  from 

constructive  interference  of  neighboring  modes  has  also  been 

3  17  3  23 

discussed  by  Tolstoy  and  Clay  *  and  Weston. 

The  work  reported  in  this  section  concerns  an  extension 
of  the  work  of  Tindle  and  Guthrie  to  a  range  dependent  medium. 
Their  techniques  are  extended  to  an  adiabatic  wave  theory  descrip¬ 
tion  of  sound  propagation  in  an  effort  to  determine  whether 


multipath  conversion  phenomena  are  contained  within  the  adiabatic 
solution  or  whether  they  are  attributable  to  mode-mode  coupling 
processes. 

What  is  meant  by  the  term  "multipath  conversion"  is 
illustrated  in  Fig.  III. 13.  Figure  III. 13  shows  how  three  rays 
initially  having  the  same  grazing  angles  with  respect  to  horizontal 
are  converted  into  three  completely  different  types  of  raypaths  by 
virtue  of  their  interactions  with  a  sloping  bottom.  With  each 
bottom  interaction  a  ray  is  converted  into  another  ray  that  is 
steeper  or  shallower  in  angle  by  twice  the  local  bottom  slope 
depending  on  whether  the  ray  is  traveling  up  or  down  slope, 
respectively.  The  multipath  conversion  process  has  a  great 

influence  on  sound  propagation  over  range  variable  bottoms  and  it 
is  important  to  know  whether  the  process  is  contained  within 
adiabatic  mode  theory  or  whether  it  was  discarded  with  the  coupling 

coefficients  and  the  assumption  of  adiabacity  of  the  medium.  If 

multipath  conversion  is  not  contained  within  adiabatic  mode  theory, 
then  the  adiabatic  approximation  would  be  of  little  use  in 

realistic  propagation  problems  involving  range  variable  boundaries. 

1.  Raypaths  Within  the  Context  of  Adiabatic  Mode  Theory 
In  analogy  to  the  work  of  Tindle  and  Guthrie  in 
range  invariant  media,  it  is  assumed  that  raypaths  linking  source 
and  receiver  correspond  to  groups  of  normal  modes  in  the  summation 
of  Eq.  (3.7)  that  interfere  constructively.  This  assumption,  which 
was  shown  to  be  correct  in  the  horizontally  stratified  case,  means 
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that  a  raypath  is  given  by  the  locus  of  range  and  depth  points  at 
which  the  phase  differences  between  a  group  of  successive  terms  in 
(3.7)  are  integer  multiples  of  2  .  Since  there  is,  in  general, 
more  than  one  raypath  connecting  source  and  receiver  locations, 
there  will  be  more  than  one  unique  set  of  r  and  z  points  that 
correspond  to  constructive  interference  of  groups  of  normal  modes 
in  (3.7). 

As  in  the  work  of  Tindle  and  Guthrie,  the  WKB 
approximations  for  the  depth  and  radial  functions  are  employed. 
These  are  given  by 


<b  (z:r)  =  const  x 
n 


/k  (z;r) 


j+i  y*  <n(z';r)dz'|  •  (3.37, 


R  (rl  =  const  x  - - - exp  \i  f  k^(r')dr'  [  .  (3.38) 

n  /TW  l  Jo  n  J 

n 

In  Eq.  (3.37)  zq  is  the  source  depth  and  the  integral  from  zq  to  z 
may  proceed  directly  or  through  the  upper  turning  point.  The 
function  <  is  the  vertical  wave  number  given  by 


(3.39) 


The  phase  of  each  term  in  the  expression  for  the  field  given  by 
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<Mr,z)  =  Vr)  4>n(z;r)  (3.40) 

n 

is 

rn(z'r)=±/*  <n(z';r)dz'  +  f  kfr'Jdr*  .  (3.41) 

zo  ° 

The  ±  sign  in  (3.41)  accounts  for  up  and  downgoing  paths  from  the 
source . 

As  stated  earlier,  a  raypath  is  assumed  to  be  given 
by  the  set  of  r  and  z  points  such  that 


Ar 

_ n 

An 


2«,tt 


(3.42) 


with  ^  being  an  integer.  The  quantity  Ar  /An  means  change  in  F 

n  n 

corresponding  to  a  change  of  mode  index  of  one  unit.  If  A/An  is 
approximated  by  9/3n,  (3.41)  in  (3.42)  gives 


/ 


z  3k 

( z ' ;  r) dz '  + 


/ 


3k 

n 

3n 


(r  *)dr '  =  2U 


From  (3.39)  it  follows  that 


3k 
_ n 

3n 


k  3n 
n 


(3.43) 


(3.44) 


hence,  Eq.  (3.43)  becomes 
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I  kn(r) 


3knlr)  fZ  dz'  , 
3„  J  Kn(z’;r) 


r  8k  (r') 
n 

8n 


dr'  =  2  J2.7T 


.  (3.45) 


o 

Equation  (3.45)  is  an  integral  equation  which  is  assumed  to 
prescribe  combinations  of  r  and  z  points  corresponding  to  rays 
emanating  from  the  source. 

To  see  whether  the  raypaths  described  by  (3.45) 
exhibit  multipath  conversion  effects,  a  test  case  is  designed.  The 
test  case  will  consider  propagation  upslope  in  an  isovelocity, 
wedge  shaped  ocean.  The  source  is  taken  to  be  at  the  surface  z*=0 
and  the  bottom  is  taken  to  have  a  linearly  decreasing  depth  given  by 


H(r)  =  Hq(1  -  Br) 


(3.46) 


In  (3.46)  B  is  a  positive  constant  for  upslope  propagation.  For  a 
rigid  sloping  bottom  and  pressure  release  surface  the  normalized 
depth  functions  are  given  by 


<t>n(z;r) 


sin 


K  Z 

n 


with 


(3.47) 


k  (r) 
n 


./ 


k2  -  kn2(r) 


(n  +  1/2)  tt 
H(r) 


(3.48) 


and  k“u/c. 
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The  test  case  chosen  here  has  an  obvious  advantage 

in  that  the  vertical  and  horizontal  wavenumbers  K  and  k  are  known 

n  n 

analytically  and  that  conventional  ray  theory  paths  are  easily 
generated  for  comparison  with  Eq.  (3.45).  From  (3.47)  and  (3.48)  it 
is  easily  seen  that 


3k  k  9k 

_ n  _  _n  _ n 

9n  "  k  9n  "  (3.49) 

n 


9k 
_ n 

9n 


IT 

H(r) 


Equations  (3.48)-(3.50)  in  (3.45)  give 


(3.50) 


r  .  1/2).  dr,  , 

0  If  frMHVrn 


2  Hit 


(3.51) 


kn(r')H  (r») 


The  integral  in  (3.51)  may  be  done  analytically  by 
changing  the  variable  of  integration  from  r'  to  H.  With  this 
procedure,  Eq.  (3.51)  becomes 


T r  .  TT  -1/Kn\ 

1  hTTY  u  -  z0)  +  Ocos  (x) 


2£tt 


Hf 


(¥)  (z  ‘  zo}  +  gF"  ton(r)  "  °n(0)1  =  2U 


(3.52) 
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a  (r)  =  cos 
IK  ' 


(3.53) 


In  terms  of  the  equivalent  ray  angle  &n(r)  defined  by 


9  (r)  =  cos 
nl 


rkn  (r)  "j 

— r —  =  o 

k  J  n 


(r)  -  tt/2 


Eq.  (3.52)  becomes 


H  ( 


V+HT[en(r)-6n(0)l  =  2^ 


(3.54) 


Equation  (3.54)  gives  the  relationship  between  r  and 
z  that  describes  the  paths  of  constructive  interference  emanating 
from  the  source  at  a  grazing  angle  of  6^(0).  The  integer  l  denotes 
the  number  of  bottom  bounces  incurred  along  the  raypath.  To  see 
whether  the  adiabatic  solution  exhibits  multipath  conversion,  the 
paths  of  constructive  interference  given  by  (3.54)  will  be  compared 
to  conventional  ray  theory  paths  for  the  same  sloping  bottom 
geometry.  When  a  ray  encounters  a  sloping  bottom  it  is  converted  to 
a  grazing  angle  that  is  steeper  by  twice  the  bottom  slope  a  as  shown 
in  Fig.  111.14.  If  Eq.  (3.54)  gives  this  same  sort  of  behavior  it 
will  be  possible  to  conclude  that  the  adiabatic  solution  contains 
within  it  the  effects  of  multipath  conversion. 

2.  Results 

In  this  section,  the  results  of  ray  tracing  via  Eq. 
(3.54)  will  be  shown.  For  simplicity  the  source  is  taken  to  be  at 


~7Vrr— 
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the  surface.  In  the  examples  to  be  presented,  the  isovelocity  wedge 
was  taken  to  be  0.5  km  deep  at  the  source,  i.e.,  Hq=0.5  km.  The 
frequency  and  sound  speed  were  taken  to  be  20  Hz  and  1.5  km/sec, 
respectively.  For  this  set  of  geoacoustic  parameters,  there  are  13 
discreet  modes  having  mode  numbers  0  through  12  present  at  the 
source. 

Ray  traces  for  two  different  launch  angles  and  three 
different  bottom  slopes  are  given.  The  two  launch  angles 
considered  are  10.81°  and  39.61°.  These  angles  correspond  to  the 
equivalent  ray  angles  of  modes  2  and  8  at  the  source.  The  three 
different  bottom  slopes  considered  are  0.5°,  2°,  and  10°.  The 
tangent  of  the  bottom  slope  angle  (a)  is  given  by 

tan  a  =  H 

o 

Figures  III. 15  through  III. 20  depict  ray  traces  performed  with  Eq. 
(3.54)  for  the  six  combinations  of  slopes  and  launch  angles.  Also 
shown  on  the  figures  are  raypaths  obtained  by  conventional  ray 
tracing  methods.  Note  the  excellent  agreement.  The  conventional 
ray  traces  were  obtained  from  a  two-dimensional  ray  trace  program 
developed  at  ARL:UT. 

In  Figs.  III. 15  through  III. 20  there  is  good 
agreement  between  the  conventional  ray  paths  and  those  calculated 
using  Eq.  (3.54),  even  for  a  bottom  slope  of  10°.  This  finding  is  a 
bit  puzzling  especially  in  light  of  the  work  presented  in  section  B 
of  this  chapter  where  it  was  found  that  the  adiabatic  approximation 
can  begin  to  be  suspect  for  bottom  slopes  as  small  as  0.5°.  Figures 
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FIGURE  111.15 

RAY  TRACE  PERFORMED  BOTH  BY  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  =  10.81°,  SLOPE  ANGLE  =  0.5° 
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CONVENTIONAL  RAY  TRACE 
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FIGURE  111.16 

RAY  TRACE  PERFORMED  BOTH  BY  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  =  39.61°,  SLOPE  ANGLE  =  0.5° 
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FIGURE  111.17 

RAY  TRACE  PERFORMED  BOTH  BY  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  =  10.81°,  SLOPE  ANGLE  =  2.0° 
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FIGURE  111.18 

RAY  TRACE  PERFORMED  BOTH  BY  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  =  39.61°,  SLOPE  ANGLE  =  2.0° 


DEPTH 


nii 


400  h 
500 


-  CONVENTIONAL  RAY  TRACE  - 

+  O  OBTAINED  FROM  ADIABATIC 

10*\  SOLUTION 

i  I  i _ l _ I _ l _ l _ I - 1 - 

0.4  0.6  1.2  1.6  2.0 


RANGE  -  km 


FIGURE  111.19 

RAY  TRACE  PERFORMED  BY  BOTH  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  -  10.81°,  SLOPE  ANGLE  =  10° 
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FIGURE  111.20 

RAY  TRACE  PERFORMED  BOTH  BY  CONVENTIONAL  METHODS  AND 
WITHIN  THE  CONTEXT  OF  ADIABATIC  MODE  THEORY 
RAY  ANGLE  =  39.61°,  SLOPE  ANGLE  =  10° 


III. 15-III.20  indicate  that  the  adiabatic  solution  describes  the 
geometric  characteristics  of  the  propagation  quite  well  for  slope 
angles  up  to  10°.  This  apparent  paradox  is  due  to  the  fact  that  the 
work  presented  in  section  B  considers  amplitude  interchange  between 
modes  while  the  work  just  described  does  not. 

Milder's  criterion  used  in  section  B  involves  the 
Bmn  which  are  proportional  to  radial  rates  of  interchange  of 
amplitude  between  modes  m  and  n.  This  amplitude  interchange  does 
not  occur  in  the  adiabatic  solution  and  hence  the  criteria  for 

adiabacity  given  in  section  B  are  based  on  amplitude  interchange 
between  modes  being  small.  The  work  just  presented  does  not 
consider  amplitudes.  It  examines  the  geometric  or  ray  properties 
of  the  adiabatic  approximation  and  has  shown  that  the  adiabatic 
solution  gives  the  proper  geometrical  properties  of  the  field  even 
for  quite  large  bottom  slope  angles.  These  findings  lead  one  to 
conclude  that  the  adiabatic  solution  describes  the  phase  of  the 

acoustic  field  quite  well  and  that  the  coupling  coefficients,  which 
are  neglected  in  the  adiabatic  approximation,  predominately 
influence  the  magnitudes  of  the  radial  distribution  of  amplitude 
among  the  normal  modes.  This  conclusion  will  be  tested  by  numerical 
calculation  and  shown  to  be  valid  in  Chapter  V. 

The  conclusions  just  reached  concerning  the 
geometrical  properties  of  sound  propagation  contained  within  the 
adiabatic  solution  were  arrived  at  through  the  use  of  an 

isovelocity,  wedge  model  for  the  propagation  medium.  It  is 


asserted  now  that  these  conclusions  are  general  because  the 
isovelocity  wedge  with  a  rigid  bottom  is  about  the  most  severe  case 
one  could  apply  the  adiabatic  approximation  to.  In  section  B  of 
this  chapter  it  was  seen  that  the  harder  the  bottom  became  in  terms 
of  characteristic  acoustic  impedance,  the  less  the  local  bottom 
slope  that  could  be  tolerated  in  the  adiabatic  approximation.  A 
rigid  bottom  is  infinitely  hard  in  terms  of  acoustic  impedance; 
hence,  the  mode  conversion  effects  in  the  test  case  considered  are 
expected  to  be  maximal  with  the  adiabatic  approximation  being  the 
most  sensitive  to  bottom  slope.  The  finding  that  the  adiabatic 
approximation  gives  the  geometric  properties  of  the  field  in  this 
case,  even  for  quite  large  bottom  slopes,  leads  one  to  conclude  that 
it  does  so  in  less  severe,  more  realistic  waveguides  also. 

Another  characteristic  of  the  isovelocity  wedge  that 
makes  it  a  severe  test  case  is  the  fact  that  all  of  the  mode  depth 
functions  have  turning  points  at  the  surface  and  on  the  bottom.  In 
a  realistic  oceanic  waveguide  the  water  sound  speed  profile  often 
traps  modes  in  the  water  column  away  from  sloping  bottom 
interfaces;  hence,  only  a  certain  portion  of  the  discrete  modes 
present  interact  with  the  bottom  and  experience  mode  conversion 
effects.  In  the  isovelocity  wedge  each  and  every  mode  present 
interacts  with  the  bottom  and  experiences  mode  conversion  effects. 


CHAPTER  IV 


COUPLED  MODE  THEORY  APPLIED  TO  PROBLEMS 
WITH  RANGE  DEPENDENT  BOUNDARIES 

In  Chapter  II  it  was  pointed  out  that  the  partial  separation 
of  the  radial  and  depth  variables  was  not  strictly  possible 
whenever  range  variable  boundaries  and/or  interfaces  were  present 
in  the  propagation  medium.  This  nonseparability  arises  from  the 
requirement  that  the  particle  velocity  normal  to  a  surface  of 
discontinuity  be  continuous  across  that  surface.  For  a  field  given 
by 


iKr.z)  =  ^  Rn(r)  <{>n(z;r)  ,  (4.1) 

n 


the  requirement  that  normal  velocity  be  continuous  entangles  the 
radial  and  depth  functions  with  the  net  result  being  that  there  is 
no  advantage  to  be  had  by  expressing  the  field  in  the  form  of  (4.1). 

As  was  pointed  out  earlier  in  Chapter  II,  a  natural 
approximation  to  be  made  when  range  dependent  boundaries  are 
present  is  to  approximate  the  normal  derivative  of  the  field  at  a 
boundary  with  the  z-derivative  and  require  this  quantity  to  be 
continuous.  For  a  waveguide  having  a  range  variable  boundary  as 
shown  in  Fig.  II. 2,  this  approximation  is  mathematically  stated  by 


z=ll(r) 


9z 


z=H(r) 


(4.2) 
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with 


(4.3) 


Equation  (4.2)  involves  an  approximation  to  the  normal  derivative 
operator  of  Eq.  (4.3)  that  is  of  zero  order  in  the  local  bottom 
slope  ti. 

In  the  derivation  of  the  mathematical  formalism  of 
coupled  mode  theory,  the  radial  functions  were  shown  to  satisfy 


R  +  -  R  +  k  2(r)  R 
m  r  m  m  m 


A  R 


mn  n 


+ 


•(4.4) 


From  Eqs.  (2.44)  and  (2.51)  it  is  apparent  that  is  of  order  H 

•  o  t 

and  A  is  of  order  H  .  This  fact  makes  the  approximation  of 
mn 

Eq.  (4.2)  suspect  in  that  terms  of  order  H  and  higher  are  dropped  in 
(4.2)  while  they  are  retained  in  (4.4).  In  other  words  the  approxi¬ 
mation  in  Eq.  (4.2),  though  necessary  to  obtain  a  workable  theory, 
is  not  a  consistent  approximation.  To  make  the  theory  consistent 
with  (4.2)  would  require  neglecting  the  right-hand  side  of  (4.4). 
This  leads  to  the  adiabatic  approximation  discussed  in  Chapter  III. 
Hence  as  things  now  stand,  the  only  consistent  application  of 
coupled  mode  theory  to  probl  ms  involving  range  variable  boundaries 
is  within  the  adiabatic  approximation. 
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Later  in  this  chapter  it  is  shown  that  if  one  ignores  the 
inconsistency  of  Eq.  (4.2)  and  uses  the  full  framework  of  coupled 
mode  theory,  the  result  is  a  field  which  does  not  conserve  energy. 
The  nonconservation  of  energy  is  caused  by  the  imposition  of 
physically  improper  boundary  conditions  on  the  field.  Therefore, 
to  apply  the  notion  of  sound  propagation  via  coupled  normal  modes  to 
problems  having  range  variable  boundaries,  a  new  theory  must  be 
constructed.  The  purpose  of  this  chapter  is  to  show  that  the 
application  of  conventional  coupled  mode  theory  to  problems  with 
range  dependent  boundaries  results  in  nonconservation  of  energy  and 
to  develop  a  corrected  coupled  mode  theory  which  can  be  applied  when 
nonstratified  boundaries  are  present. 

A.  Energy  Considerations 

In  this  sec tion  it  will  be  shown  that  applying  coupled  mode  theory 
with  the  approximation  of  Eq.  (4.1)  to  waveguides  having  range  dependent 
boundaries  leads  to  a  theory  in  which  energy  is  not  conserved.  To 
illustrate  this  point,  one  first  considers  the  time  averaged  radial 
energy  flux,  fr,  given  by  (see  Ref.  4.1,  Chapter  V) 

fr  =  \  Re(pv*)  ,  (4.5) 

where  p  is  the  acoustic  pressure  and  v^_  is  the  radial  component  of 
particle  velocity.  In  terms  of  the  velocity  potential,  4>,  p  and  v^  are 
given  by 
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P  = 


iuipiji 


(4.6) 


v 

r 


dip 

dr 


(4.7) 


Equations  (4.1),  (4.6),  and  (4.7)  combine  with  (4.5)  to  give 


f 
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()>  R*  0  + 


m  m 


(4.8) 


In  Eq.  (4.8)  the  mode  functions  <J>n  are  real  and  the  dot  symbol  denotes 
dif ferentation  with  respect  to  range. 

Equation  (4.8)  can  be  put  in  the  following  form. 
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ia)p 
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R  R* 
n  n 


R*R 
n  m 
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l  n 
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(4.9) 


To  compute  the  time  rate  of  change  of  energy  in  the  radial  direction  J^, 
Eq.  (4.9)  is  integrated  over  the  surface  of  a  cylinder  centered  about  the 
depth  axis  at  range  r=0.  (See  Fig.  IV. 1.)  This  procedure  gives 

J  =  /  2trr  f  dz 

r  J  r 

o 


(4.10) 


j  =  Jv*  (r  R*  -  R*R  \  +  y*  (r  R*  -  R*R  )b  } 

r  2  ^  n  n  n  n )  jU  \  n  m  n  m  /  nmj 

In  obtaining  (4.10),  the  orthonormality  property  of  Eq.  (2.27)  and  the 

definition  for  B  in  Eq.  (2.30)  were  employed, 
nm 


FIGURE  IZ.1 

WAVEGUIDE  GEOMETRY  FOR  DERIVATION  OF  RADIAL  ENERGY  FLUX 
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The  time  averaged  rate  of  change  of  energy  (power  flow)  given 
by  (4.10)  is  a  quantity  that  should  be  independent  of  range  in  the 
absence  of  attenuation.  The  range  independence  of  is  equivalent 
to  conservation  of  energy.  If  Eq.  (4.10)  is  differentiated  with 
respect  to  range,  the  result  is 


j  =  (r  R*  -  R*R  )  +  V  (r  R*  -  R*R  \ 

r  2  1  Lu  \  n  n  n  n  /  Z— »  \  n  m  n  m  / 

l  n  '  n,m 


nm 


r  V  (r  R*  -  R*R  ( R  R*  -  R*Rm)  B  (4.11) 

\  n  n  n  n  /  i-J  \  n  m  n  m  /  nm 

n  '  '  r>  m 


n  ,m 


+  r 


El R  R*  +  R  R*  -  R*R  -  R*R  \  B  l 
\nm  nm  nm  nm/nmf 
n,m  ) 


From  (4.4)  it  follows  that 


(r  r*  ■  R*i  )*-(rr*  ■  KK) 

\nn  nn/  r\nn  n  n  / 


E|a  (RR*-  R*R  \ 
nm  \  n  m  n  m  / 

m  v. 

PR  R*  -  R*R  )  /  .  .  vl' 

JB - HJLL  +  2  ( R  R*  -  R*Rm) 

r  V  n  m  n  m/JJ 


(4.12) 
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In  (4.17)  Che  supers :ript  A  and  S  mean  "antisymmetric  and  symmetric 

parts  of,"  respectively.  For  a  waveguide  having  one  range  variable 

interface  as  in  Fig.  II. 2,  (B  -A  )  is  obtained  from  Eq.  (2.56) 

’  nm  nm 

g 

and  B  is  given  by  (2.47). 
nm  J 

(»„.  -  \„)A  '  Apl[*ln(hm  -  *2»)  -  '  *2")L„(r) 


HP 


bL  ■  -  ~T  (1  -  0l/p2) 


z=H(r) 


(4.19) 


Equations  (4.17)-(4. 19)  indicate  that  only  when  H  is  identically  zero, 
i.e.,  horizontal  boundaries,  is  energy  conserved.  This  discovery  is  a 
rather  disconcerting  one  since  a  theory  that  fails  to  conserve  energy  is 
highly  suspect. 

At  this  juncture  it  is  reasonable  to  ask  what  happens  to 
the  energy  that  is  removed  from  the  flow  in  the  radial  direction? 
The  answer  to  this  question  is  that  the  energy  flux  removed  from  the 
radial  direction  contributes  to  a  nonzero  net  flux  of  energy  normal 
to  the  range  variable  boundary  or  boundaries.  To  illustrate  this 
point,  consider  Fig.  IV. 2.  It  will  now  be  shown  that  the  difference 
in  radial  energy  flow  through  the  surfaces  of  cylinders  of  radii  r 
and  r+dr  is  equal  to  the  net  flow  of  energy  normal  to  the  interface 
z=H(r)  between  r  and  r+dr. 


FIGURE  H.2 

WAVEGUIDE  GEOMETRY  FOR  ENERGY  FLUX  CONSIDERATIONS 
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For  an  infinitesimal  dr,  the  difference  in  the  rate  of 
energy  change  between  r  and  r+dr  is  given  by 


6J 


r  =  Jr(r)  -  Jr(r  +  dr)  =  -Jr(r)  dr 


(4.20) 


With  (4.17),  Eq.  (4.20)  becomes 


.  6  J 

*  r 


lTUiir 


„  dr["  T  C*  ‘A  )  0  +  BS  P 

2  jLj  \  nm  nm  /  ^nm  nm  nm 

Ln,m  ' 


(4.21) 


The  normal  flux  of  energy  fx  is  given  by 


fx=  \  Re[Pvj£]  =  j  Re[iwpi|»VAi|>*] 


(4.22) 


with  vx=Vx\().  With  Eqs.  (4.1)  and  (4.3),  Eq.  (4.22)  may  be 
expressed  as 


f-L  = 


(Xip 


Im 


7^{E  fa 

u.»  L 


9<J> 

A  n 
nrm  dz 


«(Wn*m  *  Wm*n) 


(4.23) 


The  net  flux  normal  to  the  interface  defined  L,  z=H(r)  is 


6f^  =  f[^  -  fj_2) 


z=H(r) 


(4.24) 
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where  the  superscript  indicates  the  layer  in  which  the  normal  flux 
is  to  be  evaluated.  Using  Eq.  (4.24)  and  the  boundary  conditions 
on  •t’n  in  Eq.  (4.24),  one  may  obtain 


(4.25) 


Since  <f> .  4> .  is 
In  1  m 


symmetr ic , 


R*R  -R  R*  may  be  replaced  with  its 
n  m  n  m  J  r 


symmetric  part  given  by 
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-  R  R* 
n  m 


R*R  -  R  R*  +  R*R 
n  m  n  m  m  n 
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R  R* 
m  n 
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nm 
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(4.26) 


With  Eqs.  (4.26)  and  (4.14),  Eq.  (4.25)  becomes 

P 


6fA  = 


/TF  n’m 


(4.27) 


'  '  ^n)} 


z=H(r) 


Since  is  antisymmetric,  Eqs.  (4.18)  and  (4.19)  combine  with 

(4.27)  to  yield 


«fX  = 


y  Bs  P  +  Q  /B  -  A  ) 
f~. j  jf-*  nm  nm  vnm\  nm  nm  / 

l/l-H  ’ 


(4.28) 


wmmmm 
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The  time  rate  of  change  of  energy  normal  to  z=H(r) 
between  r  and  r  +  dr  is  obtained  from  (4.28) 


6JX  =  6fl  6S 


(4.29) 


In  (4.29)  6S  is  the  surface  area  normal  to  z=H(r)  between  r  and  r+dr 
and  is  given  by 


6S 


2lTr  ?r  =  2 nr  dr 
cos  a 


(4.30) 


The  angle  a  in  (4.30)  is  the  local  slope  angle  of  the  surface  z=H(r) 
at  range  r,  i.e.,  H=tana.  From  (4.28)-(4.30)  <5jj_  becomes 


<SJ. 


icoiTr  dr 


n,m 


B  P  + 
nm  nm 


(B  -A  ) 
\  nm  nm  I 


(4.31) 


A  comparison  of  Eqs.  (4.21)  and  (4.31)  shows  that  they  are  the  same; 
hence,  the  change  in  radial  energy  flow  occurring  over  an  increment 
dr  is  due  to  the  net  flow  of  energy  normal  to  the  sloping  interface 
in  the  same  increment. 

The  nonzero  net  flow  of  energy  normal  to  range  variable 
boundaries  is  a  consequence  of  the  approximate  boundary  conditions 
imposed  on  the  field.  The  mathematical  formalism  of  coupled  mode 
theory  determines  a  field  satisfying  precisely  the  boundary 
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conditions  it  was  requested  to  (continuity  of  pressure  and  z- 
component  of  particle  velocity).  Should  such  a  medium  exist  that 
required  continuity  of  pressure  and  z-component  of  particle 
velocity,  conventional  coupled  mode  theory  would  produce  the 
correct  field  and  radial  energy  flow  would  not  be  range 
independent.  However,  when  this  field  is  used  to  approximate  a 
field  which  is  supposed  to  satisfy  the  physically  proper  boundary 
conditions  (continuity  of  pressure  and  normal  particle  velocity)  the 
the  nonconservation  of  energy  byproduct  is  unacceptable. 

The  ultimate  importance  of  the  nonconservation  of  energy 
effect  in  realistic  applications  of  coupled  mode  theory  depends  on 
its  importance  compared  to  any  attenuation  effects  expected  to  be 
present  in  the  problem.  An  estimate  of  the  relative  importance  of 
energy  and  attenuation  effects  can  be  made  as  follows.  Equations 
(4.10)  and  (4.17)  may  be  expressed  as 


(4.32) 


(4.33) 


with 
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i  TToor 
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(4.34) 
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nm 


(4.35) 


In  Eqs.  (4.34)  and  (4.35)  P  and  Q  are  given  by  (4.14)  and 

nm  nm 

(4.15).  The  may  be  thought  of  as  the  time  rate  of  change  of 

energy  in  the  radial  direction  associated  with  mode  number  n. 

The  ratio  of  Eq.  (4.35)  to  (4.34)  to  first  order  in  the 
bottom  slope,  H,  is  given  by 
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nm 
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nm 
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nn 
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0  (H  ) 


(4.36) 


The  quantity  $n  may  be  thought  of  as  an  effective  power  attenuation 
coefficient  describing  the  reduction  in  radial  energy  flow  due  to 
the  nonconservation  of  energy  effects.  From  Eq.  (4.36)  it  follows 
that 

/n)  (r)  =  J^n)  (rQ)  exp Sf  iyr')dr'J  (4.37) 

^  r0 

with  j£n^(rQ)  being  the  time  rate  of  energy  flow  associated  with 
mode  n  at  some  reference  range  rQ. 

The  attenuation  of  the  power  associated  with  mode  n  may 
be  estimated  from  the  results  of  Chapter  II.  From  Eqs.  (2.69)  and 
(2.70),  the  imaginary  part  of  the  eigenvalue  resulting  from 


attenuation  in  the  medium  is  seen  to  be 
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i6n(r) 

In>[kn(r)]  »  ‘  (4'38) 

The  attenuation  to  the  radial  functions  R^  caused  by  the  imaginary 
part  of  can  be  approximately  included  by  replacing  R^  with 

f  rT  6n(r,) 

Rn  exp  |  "  jf  2k-(r’7  dr’j  •  (4-39) 

Since  the  power  associated  with  mode  n  is  quadra tical ly  related  to  the 
amplitude  R^,  the  power  attenuation  due  to  absorption  in  the  medium  is 
expressed  approximately  by 

r  rT  <Vr'-)dr'  1 

exp fir  •  <4-40> 

^  0 

A  comparison  of  Eqs.  (4.40)  and  (4.37)  indicates  that  one  can  expect 
the  nonconservation  of  energy  effect  to  be  negligible  compared  to 
normally  occurring  attenuation  when 

r  rr  6  (r')dr' 

J  3  (r')dr'  «  J  -  .  (4.40 

*>  r  n 

o  o 

Since  Eq.  (4.36)  is  good  to  0(H)  only,  one  is  entitled  to  use  the 

adiabatic  radial  solutions  in  the  evaluation  of  6  .  Therefore  a 

n 

determination  as  to  the  importance  of  nonconservation  of  energy 
effects  based  on  Eq.  (4.41)  may  be  made  using  adiabatic  mode  theory 
quantitites  which  are  relatively  easy  to  obtain. 
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The  nonconservation  of  energy  phenomenon  has  not  arisen 
before  in  the  relatively  few  applications  of  coupled  mode  theory  to 
acoustic  propagation  problems  because  most  of  the  applications  have 
involved  horizontally  stratified  boundaries  and  the  rest  failed  to 
examine  or  consider  energy  conservation.  In  Refs.  4. 2-4. 4, 
the  coupled  mode  theory  formalism  was  applied  to  underwater 
acoustic  propagation  problems  involving  range  dependent  sound  speed 
with  either  no  boundaries  or  horizontally  stratified  boundaries. 
In  Refs.  4. 5-4. 7,  the  existence  of  nonhorizontal  boundaries 
was  allowed  for,  but  energy  conservation  questions  were  not  considered. 

B.  First  Order  Correction  to  Coupled  Mode  Theory 

In  this  section  a  correction  to  the  original  coupled  mode 
theory  of  Pierce  and  Milder  is  proposed  and  derived.  It  will  be 
shown  that,  with  rather  modest  change,  coupled  mode  theory  may  be 
applied  to  problems  with  range  variable  boundaries  with  the  field 
so  obtained  satisfying  the  proper  boundary  conditions  and 
conserving  energy,  both  to  first  order  in  the  boundary  slopes  H. 

To  derive  the  appropriate  corrections  to  coupled  mode 
theory,  Eq.  (4.1)  is  replaced  with 

ip(r,z)  =^[RnO)  +  <5Rn(r)]  [4>n(z;r)  +  6<)>n(z;r)]  .  (4.42) 

n 

In  (4.42)  the  6  R  and  66  are  correction  terms  of  first  order  in 
n  n 

are  included  so  that  ip  can  be  made  to  satisfy  the  physically  correct 
boundary  conditions  to  first  order.  The  $  and  Rn  of  Eq.  (4.42) 
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are  the  same  as  those  described  in  Chapter  II.  The  6  Rn  is 

taken  to  be  solely  a  function  r  with  6R  and  all  its  range 

n 

derivatives  being  of  order  H.  The  64>n  is  the  corrective  term  added 
to  4>n.  Like  <t>n(z;r),  6<J>n(z;r)  is  predominately  a  function  of  z  and 
is  weakly  dependent  on  range. 

At  this  point  it  is  convenient  to  list  various  quantities 
that  will  arise  as  to  their  order  in  powers  of  H. 
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In  (A. 43)  the  dot  and  prime  symbols  stand  for  differentiation  with 

respect  to  r  and  z,  respectively.  The  corrections  to  coupled  mode 

theory  being  proposed  in  this  section  are  of  first  order.  Therefore 

•  2  . 

any  terms  of  order  H  and  higher  and  combinations  of  terms  from 

•  2 

(4.43)  that  are  order  H  and  higher  are  not  carried  along. 

1 .  Boundary  Conditions  on  6^ 

The  boundary  conditions  to  be  satisfied  by  the  6$ 
may  be  obtained  by  considering  the  boundary  conditions  on  The 
requirement  that  pressure  be  continuous  everywhere  means 
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that  P*P  must  be  continuous  at  boundaries  and  interfaces.  The  6 

n 

are  defined  in  Chapter  II  such  that  p4>  is  a  continuous  function  of 

n 

depth.  Therefore,  from  (4.42),  if  P'11  is  to  be  a  continuous 

function  of  depth,  p<5</>  must  be  required  to  be  continuous.  In 

n 

addition,  since  the  field  must  vanish  at  0  and  °°,  6(f)  must  vanish  at 

n 

0  and  00  also. 

The  particle  velocity  normal  to  an  interface  or 
boundary  is  also  required  to  be  continuous.  From  Eq.  (4.3),  the 
particle  velocity  normal  to  an  interface  defined  by  z=H(r)  is  to 
order  H. 


V  ■  H  - » k*  ♦  °W2> 


(4.44) 


Using  Eq.  (4.42)  in  (4.44)  one  obtains  to  first  order 


E(r  +  6R  )  0'  +  R6<P'  -  HR  *  1  +0(H2) 
\n  n/n  nn  nnl 

n  '  J 


(4.45) 


where  (4.43)  has  been  employed.  Equation  (4.45)  is  to  be  continuous 

at  surfaces  of  discontinuity.  Since  <J>  is  defined  such  that  ’  is 

n  n 

continuous  everywhere, 


R  6<P'  -  HR  <(> 
n  n  nTn 


(4.46) 


must  be  required  to  be  continuous.  For  a  nonhorizontal  boundary 
defined  by  z=H(r)  separating  two  fluids  of  densities  and  p2  (see 
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Fig.  II. 2)  the  discontinuity  in  <5<J>'  at  z=H(r)  obtained  from  (A. 46) 

n 

is 


z=H(r) 


R  . 

if  H*ln  (1-Pl/P2) 


z=II(r) 


( 4.47  ) 


In  (4.47)  the  subscript  preceding  the  mode  index  denotes  the  layer 
in  which  the  corresponding  function  is  to  be  evaluated. 


2.  The  Differentia?  Equation  for 

An  equation  for  the  SR^  may  be  obtained  by  using  Eq. 
(4.42)  in  the  differential  equation  satisfied  by  \p . 


V2ij> 


+  k  (r,z)ip=0 


(4.48) 


Equation  (4.42)  in  (4.48)  yields  to 
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In  deriving  Eq.  (4.49)  the  fact  that  £  Rn4>n  satisfies  (4.48)  and  the 
equation  for  <p ,  (2.24),  was  used.  If  (4.49)  is  multiplied  by  p<J>m 


IWPliWyWi 


and  integrated  over  depth  and  if  the  orthonormality  property  of  Eq. 
(2.27)  is  employed,  (4.49)  becomes 
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(4.51) 


The  second  term  on  the  right-hand  side  of  (4.50)  may 
be  integrated  twice  by  parts  with  respect  to  z.  If  one  recognizes 
that  the  integral  from  0  to  00  has  the  meaning  given  in  Eq.  (2.34), 
the  following  expression  is  obtained. 


■  (,2*2mS»2n]H(r) 

[°l*i„6*l„  -  p2^2m^2n]  .  , 
1  JH(r) 


(4.52) 
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+ 


6R 


m 


where 


R  I 
m  mn 


2BS  R 
mn 


(A. 56) 


mn 


(4.57) 


from  Eq.  (2.47)  has  been  used. 

Equation  ^4.56),  though  derived  with  a  particular 
waveguide  structure  in  mind  (e.g.,  Fig.  II. 2),  is  completely 
general  in  form.  If  one  desired  to  include  more  than  one  non¬ 
horizontal  interface,  Eq.  (4.53)  would  have  discontinuity  terms  at 
each  interface  as  would  (4.57),  the  net  result  being  that  (4.56) 
would  remain  identical  in  form.  The  derivations  presented  here  are 
based  on  one  sloping  interface  defined  by  z=H(r)  solely  to 
illustrate  the  first  order  correction  technique  without  unduly 
complicating  the  already  involved  analysis.  Before  further 
simplification  of  Eq.  (4.56)  can  be  carried  out,  the  differential 

equation  for  6$  must  be  obtained  so  that  I  may  be  evaluated. 
n  n  mn  7 

3.  The  Equation  for  6<t>n 

To  obtain  an  equation  for  <5<t>n»  °ne  assumes  that  the 
combination  <^n=^1+^<^n  satisfies  an  equation  very  similar  to  the 
one  satisfied  by  e.g., 


9  ^ 
— j  +  k  (r , z) 

Laz 


cn2(r) 


<J>  = 

n 


0 


(4.58) 
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In  (4.58)  *7  is  an  altered  eigenvalue  which  will  be  shown  to  differ 
n 

from  k  by  terms  of  order  H.  Since  0'  is  continuous  everywhere,  0  ' 
n  n  }  ’  n 

will  satisfy  the  same  boundary  conditions  as  <$<}>  f  i.e., 


'in  -  *2n  '  ¥\.(«  -  V0?) 


z=H(r) 


(4.59) 


p.  4>,  = 

1  In  2  2n 


z=H(r) 


The  assumption  of  the  form  of  (4.58)  is  based  on  the  requirement 

that  ^  approach  as  H  goes  to  zero.  This  must  be  true  because 

the  conventional  coupled  mode  theory  is  valid  for  H=0.  From  Eq. 

(4.59)  it  is  apparent  that  the  boundary condi t ions  on  ♦  become  the 

same  as  those  on  as  H  goes  to  zero.  With  Eqs.  (4.58)  and  (4.59) 

one  i3  assured  that  the  corrected  theory  reduces  to  conventional 

coupled  mode  theory  in  the  limit  of  horizontal  boundaries. 

The  equation  satisfied  by  <5<|>n  is  obtained  by  substituting 

the  equation  for  <P  »  which  is 
n 


+  k  (r,z) 


k  2(r) 
n 


4)^  C  z ;  r)  =  0 


(4.60) 


into  (4.58).  This  procedure  yields 


(^'k2-°n2)6*n-K2-kn2)*n 


(4.61) 


•  .  .  .  • 
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X 


X 


/ 


It  is  now  possible  to  derive  an  expression  for  I  .  From  (4.60)  and 


(4.61)  the  following  expression  may  be  obtained. 


<t  6<b"  -  6d>  <b"  =  (a  -  k  “)<J>  6cf>  +  (o  '  -  k  "W  <J 
Tn  Ynym  \n  m/mn  \n  n'n 


(4.62) 


If  (4.62)  is  multiplied  by  P  and  integrated,  the  following  is 
obtained. 


[VXin  -  p2»X2"]H(r)  ■  [xx-  ■ 

=  (°n  ‘  *Sn  )  ('inn  +  ^n,m) 


r) 


(4. S3) 


If"  the  boundary  conditions  satisfied  by  and  are  used,  Eq. 

(4.63)  becomes 


p  d>  ["od)'  —  64>  ’  ~1  =  (a  “  -  k  +6 

Mlvlm[  yln  ^^n-JH(r)  '  n  m  / 1  mn  n,mj 


Equations  (4.47)  and  (4.57)  combine  with  (4.64)  to  yield 


(4.64) 


2R 

\l  +  6  l(°  2  '  k  2)  =  '  -n 

[  mn  n,mj\  n  m  /  R 


2Rn  c 
n  BS 

mn 
n 


(4.65) 


Two  important  pieces  of  information  may  be  obtained 
from  Eq.  (4.65).  For  m=n,  (4.65)  gives 


■'■“ST- 
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2  .  2 
a  -  k 
n  n 


2  .  2 
a  -  k 
n  n 


2R  , 

B  _ * — - 

Rn  nn  (1+'nn) 

2R 

R~  Bnn  +  0(H2> 
n 


(4.66) 


The  second  of  Eqs.  (4.66)  is  obtained  because  Bnn  and  Inn  are 
order  H.  Equations  (4.66)  state  that  the  difference  between  the 
eigenvalue  k  and  the  corrected  eigenvalue  is  first  order  in  H, 
verifying  that  0  approaches  k  as  H  goes  to  zero.  For  m^n  (4.65) 


gives 


-7  T  v  2R 

(a  2  -  k  2  )  1  =  -  -p11 

\  n  m  /  mn  R 

'  n 


mn 


(4.67) 


2  2 

Since  the  difference  between  0^  and  k^  is  first  order  as  is  Imn> 
Eq.  (4.67)  becomes 


(k  2  -  k  2  )  I 
\  n  m  /  mn 


2  R  ~ 

BS  +  0(H2) 
R  mn 
n 


(4.6P) 


for  m^n. 


4.  Final  Expressions  for  the  First  Order  Correct  Field 

Now  that  the  equations  for  6  R  ,  <S(J)  ,  and  I  are 
1  m  m  mn 

known,  the  component  parts  of  the  theory  may  be  assembled  in  a  form 
that  is  computationally  convenient.  Equation  (4.42)  is  best 
expressed  in  the  form 


4>(r,z)  =  £  xn(r)  fon(z;r)  +  6<J>n(z;r)] 


(4.69) 
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with 

,  Xn(r)  =  Rn(r)  +  6Rn(r)  ‘  (4.70) 

The  equation  for  is  obtained  by  adding  the  equations  for  Rr 

and  6R  given  in  (4.56)  and  (4.4).  This  gives 
n 

n  l 

(4.71) 

With  Eq.  (4.68),  (4.71)  becomes 


-  +  km2xm  =  - 
r  mm 


In  Eq.  (4.72),  if  a  term,  Rn/r,  which  is  negligible  compared  to  2Rr 
for  r  greater  than  a  few  wavelengths,  is  added  to  the  second  and 
third  terms  in  (4.72),  an  equation  comparable  in  form  to  (4.4)  is 
obtained.  Hence,  Eq.  (4.72)  becomes 


X  + 
m 


X 

m 

r 


+ 


m  m 


VB  ( ♦  2R  ) 
Cm*  mn\  r  n  / 
n*m 


(4.73) 
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In  deriving  (4.73)  the  second  order  terms  in  (4.4)  involving  A 

ran 

were  neglected.  Since  R  =X  -<$R  differs  from  X  by  terms  of  order 

n  n  n  n 

H,  one  is  entitled  to  write  (4.73)  as 


X  + 
m 


k  2X 
m  m 


VB  ( ♦  2X  ) 
mn\  r  n  / 

n=m 


(4.74) 


Equation  (4.74)  differs  from  (4.4)  in  first  order  by  the  exclusion 

of  the  m=n  term  from  the  summation.  When  H=0,  Ea .  (4.57)  shows  that 

B  is  identically  zero;  hence,  the  equation  for  X  in  the  limit 
nn  7  ’  1  n 

becomes  the  same  as  the  one  for  R  as  required. 

n 

The  equation  for  given  by  f 4 . 6 1 )  along  with 

(4.66)  becomes 


. 


.  2  2 
~-  +  k  -  a 
2  n 


2R 

6<f>  =  -  —  B  4> 

rn  R  nn  n 

n 


(4.75) 


As  stated  before,  the  desired  solution  to  Eq.  (4.75)  is  the  one 

which  vanishes  as  H  goes  to  zero.  The  homogeneous  solution  to 

(4.75)  is  not  the  proper  one  because  it  approaches  <J>  as  H  goes  to 

zero.  The  particular  solution  of  (4.75)  is  the  proper  solution  that 

vanishes  in  the  limit  since  the  right-hand  side  of  (4.75)  vanishes 

for  H=0.  Equation  (4.75)  is  not  particularly  useful  for 

computing  6$  since  it  involves  the  solution  of  an  inhomogeneous 

eigenvalue  problem.  An  approximation  to  (4.75)  involves 

2  2 

replacing  0^  with  k^  .  One  can  make  this  replacement  because  6$^ 

7-  2  . 

is  first  order  as  is  the  difference  between  0  '-k  .  With  this 

n  n 

substitution  (4.75)  becomes 
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9  v2  , 

— j  +  k  -  k 

9z 


2R 

n  ]5*n  =  -  V1 
J  n 


B  d> 
nn  n 


(4.76) 


Equation  (4.76)  is  much  easier  to  solve  than  (4.75)  because  it  does 

not  involve  computing  eigenvalues,  cr^  .  Once  the  <J>n  and  are 

known,  a  straightforward  integration  of  (4.76)  produces  the  6<j>  . 
2 

The  are  not  required  to  compute  the  corrected  field;  hence,  a 

savings  in  computational  effort  is  to  be  gained  by  using  (4.76)  over 


(4.75). 


To  show  that  <$4>^  computed  from  (4.76)  still 


satisfies  (4.47),  the  following  equation  formed  from  (4.76)  and 
(4.60)  is  considered. 


[<J>  64> ’  -  6<(>  <J»'] 
9z  1  m  n  n  m 


— —  B  4)  <f>  + 

R  nn  n  m 
n 


2  -  k  2 
i  m 


/  m  n 


•(4.77) 


If  (4.77)  is  multiplied  by  the  density  and  integrated  over  depth, 
(4.77)  becomes 


^n]2=H(r) 


-5- -  B  <5  + 

R  nn  n,m 
n 


(k  2  -  k  2  )  I  .  (4.78) 
\  n  m  /  mn 


For  m=n,  (4.78)  gives 


’“♦In  -  “♦in]2s 


z=H(r) 


2Rn  Bnn 

Rn  P^lm  z=H(r) 


(4.79) 
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With  Eq.  (4.57),  (4.79)  is  seen  to  be  the  same  as  (4.47). 

To  summarize,  a  correction  to  coupled  mode  theory 

that  is  of  first  order  in  the  slopes  of  nonhorizontal  boundaries 

present  in  the  waveguide  has  been  derived.  The  corrected  field  is 

given  by  Eq.  (4.69).  The  xn»<J>n>  and  60 n  satisfy  Eqs.  (4.74), 

(4.60),  and  (4.75),  respectively.  In  the  limit  of  H=0,  Xr  reduces  to 

and  60  goes  to  zero  such  that  the  conventional  coupled  mode 

theory  formalism  is  regained. 

The  main  difference  between  the  conventional  theory 

and  the  corrected  theory  just  derived  is  the  additional 

quantity  60  that  must  be  computed.  The  correction  to  the  radial 
n 

function  6Rr  has  been  incorporated  into  the  equation  for  X^  and 
need  not  be  computed  independently.  Depending  on  the  type  of  appli¬ 
cation,  the  added  computational  difficulty  imposed  by  the  64>n  may 
be  negligible  compared  to  the  other  calculations  that  must  be 
performed,  e.g.,  solving  Eq.  (4.74).  The  60n  are  not  needed  until 

the  field  is  to  be  constructed  from  the  X  ,  <b  ,  and  64)  .  This  is 

n’  Tn’  Tn 

not  true  of  the  4>  which  must  be  known  as  a  function  of  range  so 
n 

that  the  coupling  coefficients  may  be  determined.  Therefore,  in 
applications  of  the  theory  in  which  one  is  interested  in  the  field 
as  a  function  of  depth  at  one  particular  range,  only  the  64>n 
corresponding  to  that  range  need  be  calculated  with  a  minimal 
increase  in  computational  effort.  If  the  field  as  a  function  of 
depth  and  range  is  desired,  then  the  computation  of  the  60 n  becomes 
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1 35 


where  F  =,/r  R  and 
m  v  m 

a  B  -  B 

„A  nm  mn  ..  DS 

D  - - ~ - =  D  "  P 

nm  2  nm  nm 


(4.83) 


The  time  averaged  radial  energy  flux,  f^_,  may  be 
obtained  by  substituting  the  field  given  by  (4.80)  into  Eqs.  (4.5)  - 
(4.7). 


If  terms  higher  than  first  order  are  dropped,  Eq.  (4.84)  becomes 


f  =  -  “  Im 
r  2r 


+  G*G 
m 


(q  g*  -  1  (d>  0  +  ct>  6ct>  +  <l>  •$<]>*  \ 

1  /  -«  (nm  2r  I  \  ®  n  m  n  n  m/ 

<t>  <j>  l  • 

n  n  m  r 


(4.85) 


As  before,  the  time  rate  of  change  of  energy  in  the  radial  direction 
is  obtained  by  integrating  (4.85)  over  a  cylindrical  surface  of 


radius  r. 
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J  =  -  coirlm 
r 


{z(v;  -  (v 


+i*  +  i  ) 

nm  mn/ 


+  G*G 
m 


.  j 

n  nmj 


(4.86) 


In  the  following  developments  certain  zero  order 
approximations  for  and  are  employed.  To  zero  order  in  H,  the 
equations  for  R  and  G  are  given  by  the  left-hand  sides  of  Eqs. 
(4.4)  and  (4.82). 


R 

R  +  —  +  k  2R  =  0  +  0(H) 
nr  n  n 


(4.87) 


G  + 
n 


0  +  0(H) 


(4.88) 


The  solutions  of  Eqs.  (4.87)  and  (4.88)  for  certain  purposes  will  be 
taken  to  be  the  WKB  solutions  (see  Refs.  4.8-4.10)  given  by 


R  (r)  =  — — —  exp 

n  /FT 

n 


ji  f  kn(r')dr'| 


0(H) 


(4.89) 
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(4. <»n) 


Equations  (4.89)  and  (4.90)  are  useful  for  replacing  radial 
derivatives  of  and  with  the  following  expressions 


R  =  ( ik  -  ~  )  R  =  ik  R  in  farfield 
n  \  n  2r/n  nn 


(4.91) 


G  =  ik  G 
n  n  n 


(4.92) 


Equation  (4.68)  along  with  (4.91)  reveals  that 


Imn+l*n  Pure  imaginary  given  by 


2i  B 

it*  mn 

I  ■+  I*  =  -  - - - 

mn  nm  k  -k 

n  m 


(4. *3) 


With  the  knowledge  that  Imn+I*ni  is  P,ire  imaginary,  Eq.  (4.86)  can  be 
written  as' 

-  tfn)  *  E  (W  *  *  '!.) 

In  nJm 

L.  v  1  (4.94) 

(G  G*  -  G*G  )B  } 

\  n  m  n  m/  nm  I 
n  ,m  J 


In  (4.94)  the  terms  G*G  /2r  were  neglected  in  comparison  to  G  G*. 

m  n  n  m 

Since  the  second  and  third  terms  of  Eq.  (4.94)  have  multiplying 
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factors  of  order  H,  the  zero  order  approximation  for  of  (4.92) 
may  be  used. 


j  =  G*  -  G*G  )  +  Y"  (G  G*  -  G*G  ) 

r  2  1  i—i  \  n  n  n  n  /  i—t  \  n  m  n  m  / 

In  m ,  n 


(4.95) 


ik  (\*  +  I 

m\  nm  mn 


To  show  that  Eq.  (4.95)  is  constant  in  range  to 
first  order,  Eq.  (4.95)  is  differentiated  with  respect  to  range  to 


~  T~  (G  G*  -  G*G  ) 

2  1 t—j  dr  \  n  n  n  n / 


+  V]  -  k  )(G  G*  +  G*G  \[b  -  ik  / 1  +  I*  1 

\  n  m/  \  n  m  n  m/[  nm  m\  mn  nm  /Jj 


(4.96) 


In  obtaining  (4.96),  Eq.  (4.92)  was  employed  and  terms  of  order  H 
were  not  differentiated  since  they  give  rise  to  second  order 
quantities.  If  Eq.  (4.82)  is  multiplied  by  G*  and  subtracted  from 
its  complex  conjugate,  the  first  term  in  (4.96)  is  seen  to  be 


4-  (G  G*  -  G*G  )  =  -  V*  2BA  (f*G  -  F  G*\ 
dr  \  n  n  n  n  J  nm  \  m  n  m  n  / 

m 

Elk  2  -  k  2  )fF*G  +  F  G*  )  I 
\m  n/\mn  mn/i 


(4.97) 


ar^n*»' 
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The  zero  order  equation  for  F^  is  the  same  as  the  one  for  (see 

Eqs •  (2.31)  and  (4.88));  therefore,  since  each  term  in  (4.97)  has  a 

first  order  multiplicative  factor,  F  may  be  replaced  with  G  . 

r  m  m 


A  (c  G*  -  G*G  )  =  -  Y'  2BA  (g*G  -  G  G*) 
dr  \nn  n  n/.  nm\  m  n  m  n ! 

ITT 

.  y  (k  2  -  k  2Vg*g  +  g  g*) 

Lu  \  m  n  /V  m  n  m  n  / 


(4.98) 


nm 


Equation  (4.93)  in  (4.96)  employing  (4.9?)  gives 


2  i  k  BA 
m  nm 


_  lTTo)  I V  / g  G*  +  G*G  ^  Ti  (k  -  k  )b  +  2i 
r  2  1  \  n  m  n  m  / 1  V  n  m/  nm 

In, m  L 

+  km(kn  '  km)('nm  +  ‘nm)  ‘  (km"  '  kn~  )  ‘nm]}  ' 


(4.99) 


Since  G  G*+G*G  is  symmetric  with  respect  to 
n  m  n  m 

interchange  of  m  and  n,  the  first  two  terms  in  (4.99)  become 


V  Yg  G*  +  G*G  )  i(k 
L-j  \  n  m  n  m  /  '  n 

n,  m 


k  )  B  +  2 i k  BA 
m '  nm  m  nm 


L(g  G*  +  G*G  'i  i(k  +  k  )  BA  =0 
\nm  nm/\n  m '  nm 

n  ,m 


Hence,  (4.99)  becomes 
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j  =  (G  G*  +  G*G  )\k  (k  -  k  Wl*  +  I  ) 

r  2  Z— i  \  n  m  n  m/  m\  n  m/V  nm  mn/ 

Ln  ,m  L 


-  (k  2  -  k  ^  I  ll 

\  m  n  /  nm  f 

With  Eqs.  (4.68),  (4.93),  and  (4.91),  Eq.  (4.100)  becomes 


(4.100) 


ITUjJ 
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i  V*  (g  G*  t  G*G  ) 
]  t—t  \  n  m  n  m/ 
Ln,m 


k  (k  -  k  ) 
m  n  m 


k  -  k 
n  m 


-  k. 


f 


“mr 


=-2_L  B  1^  =  0  (4.1 0 U. 


Therefore 


Jr  -  0  +  0fH2j  '  (4.102) 

Equation  (4.102)  is  the  desired  result  showing  that  the  time  rate  of 
change  of  energy  in  the  radial  direction  has  been  made  constant  to 
first  order  by  inclusion  of  the  corrections  to  the  theory  derived  in 
this  chapter.  In  the  next  chapter  numerical  calculations  are 
presented  comparing  conventional  coupled  mode  theory  and  the 
corrected  theory.  These  calculations  will  be  seen  to  verify  the 
mathematical  analysis  presented  in  this  chapter. 


CHAPTER  V 


NUMERICAL  CALCULATIONS 

In  this  chapter  some  coupled  mode  theory  numerical 

calculations  are  presented.  The  purpose  for  making  these  calcu¬ 
lations  is  to  illustrate  how  coupled  mode  theory  is  implemented  and 
to  verify  the  assertions  and  corrections  proposed  in  Chapter  IV. 

The  material  in  this  chapter  will  be  presented  in  the 

following  format.  In  section  A  the  waveguide  model  assumed  for  the 

calculations  is  described.  In  section  B  the  numerical  procedures 
that  are  employed  to  calculate  the  field  are  discussed.  This 

discussion  includes  the  numerical  technique  used  to  solve  the 
second  order,  coupled  radial  equations,  initial  conditions, 
convergence  criteria,  and  computer  run-time  statistics.  In  section 
C  the  numerical  calculations  are  presented.  Numerical  calculations 
of  power,  radial  functions,  depth  functions,  and  the  acoustic  field 
are  presented. 

A.  Model  Description 

The  waveguide  geometry  assumed  for  the  calculations  to  be 

presented  in  this  chapter  is  shown  in  Fig.  V.l.  Figure  V.l  is  an 

isovelocity,  wedge  shaped  model  of  an  oceanic  waveguide.  The  sound 

speed  in  the  water  is  taken  to  be  isovelocity  and  the  density 

constant  in  depth  and  range.  For  the  calculations  to  be  presented 

3 

here  the  sound  speed  and  density  were  taken  to  be  1.5  km/sec  and  10 
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kg/m  ,  respectively.  The  source  is  assumed  to  he  harmonic  in  time 
with  a  frequency  of  20  Hz  and  located  at  a  depth  of  zo=0.2  km. 

The  surface  of  the  waveguide  is  taken  to  be  pressure  release  and  the 
bottom  is  assumed  to  be  rigid.  The  depth  of  the  waveguide  is  range 
dependent  and  is  given  by 


H(  r) 


0  <  r  <  r 

s 

tan  a (r-r  )  r  >  rs 


(5.1) 


with  r  =1.0  km. 
s 

The  water  depth  in  the  source  region  is  taken  to  be  Hq=0.5  km.  The  flat 
portion  of  the  bottom  in  the  source  region  is  included  so  that  the 
initial  conditions  on  the  field  in  the  range  variable  region 
(r>rg)  might  be  obtained  from  the  analytically  known  field  in  the  source 
region  (r£rg).  In  Table  II  the  appropriate  geoacoustic  parameters  are 
collected  in  one  place. 

In  the  following  sections  of  this  chapter  the  field  ^  given  by 


n 


/r 


[<}>n(z;r)  +  6<J>n ( z ; r )  ] 


(5.21 


is  calculated.  For  the  waveguide  of  Fig.  V.l  the  ^  functions  are  known 
analytically.  The  depth  functions  satisfy 


(z;r) 


0 


(5.3) 
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TABLE  II 

GEOACOUSTIC  PARAMETERS  OF  WAVEGUIDE  DEPICTED  IN  FIGURE  V.l 

P  *  10^  kg/m^ 


c  =  1.5  km/sec 


H  =  0.5  km 

o 


z  =  0.2  km 

o 


f  =  20  Hz 


r  *  1.0  ko 

3 
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For  a  constant  k=u)/c,the  depth  functions  are  given  by 


h  =  a  sin  k  z  +  h  cos  k  z 
n  n  n  n 


iS. 4) 


k  =  fv2  -  k  2 

n  v  n 


(5.5) 


For  a  pressure  release  surface  and  rigid  bottom,  in  (5.4)  must  be 

3«n 

set  to  zero.  The  rigid  bottom  requires  that 
bottom;  hence, 


vanish  at  the 


■  r(n+l/2)i:z  1 
|,nU;r)  -  an  - j 


with 


K  fr)  =  ("*1/2)1 
Kn  r  H(r) 


(5.6) 
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Since  there  are  no  surfaces  of  discontinuity  between  the  surface 
and  bottom  of  the  waveguide,  the  depth  functions  may  be  normalized 
such  that 


/ 


H(r) 


WZ  =  1 


(5.8) 


Hence,  the  orthonormal  d>  are  given  by 

n 


Vz; 


,  rr~  .  r(n^i/2)Tizi 
,r)  "  V  H(r)  l  H(r)  J 


(5.9) 


With  the  field  given  by  (5.2)  an  assignment  of  units 

among  the  radial  and  depth  functions  must  be  made.  It  is  convenient 

to  take  the  normalization  integral  of  Eq.  (5.8)  as  unitless.  With 

this  convention,  the  units  of  the  (j>n  and  64>n  are  inverse  square 

root  of  length  or  l/*^.  The  units  of  the  velocity  potential,  'P,  are 

2 

length  square  over  time  or  L  /T;  hence,  with  the  unit  assigned  to  $ 

and  6<f>n,  Gn(r)  must  have  units  of  length  cubed  over  time  or  L  /T. 

This  unit's  convention  is  followed  throughout  this  chapter. 

The  calculations  described  in  the  next  section  require  a 

knowledge  of  the  coupling  coefficients  B  .  Since  the  $  are  known 

mn  n 

analytically  in  the  case  considered  here,  the  coupling  coefficients 
may  also  be  determined  analytically  using  the  techniques  of  Chapter 
II.  They  are  seen  to  be 


B 

mn 


(r) 


(n+l/2)2H 

H(r) 


f  -1:U 


m+n 


(m+n+1)  (m-n) 


(5.10) 
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The  antisymmetric  and  symmetric  parts  of  B  are,  respectively, 


bA  =  V  '  JSL  =  (-!  )m+n  || 

mn  2  H 


(n  +  1/2)  2  +  (m+1/2) 2 
(m+n+1 ) (m-n) 


(5-13) 


and 


mn 


B  +  B 
mn  nm 

2 


(-D 


m+n+1  H 
H 


(5.12) 


In  Eqs.  (5- 10)  —  (5. 12)  the  density  is  missing  because  of  the 
slightly  altered  orthonormality  condition  of  Eq.  (5.8). 

B.  Numerical  Procedures 

In  this  section  the  numerical  procedures  that  are  used  to 

compute  the  acoustic  field  in  the  waveguide  of  Fig.  V.l  using 

coupled  mode  theory  are  discussed.  Since  the  <t>n  and  are  known 

analytically  in  this  case,  two  steps  in  the  numerical  computations 

are  bypassed.  With  the  field  given  by  (5.2)  two  quantities  remain 

to  be  calculated,  the  G  and  6<j>  .  It  turns  out  that  an  analytic 

solution  for  6<J>n  is  possible  in  the  present  case;  therefore, 

discussion  of  the  computation  of  64>  is  postponed  until  section  C. 

n 

As  indicated  in  previous  chapters,  the  real  challenge  in  coupled 

mode  theory  concerns  the  solution  of  the  coupled  radial  equations. 

This  section  focuses  on  that  aspect  of  the  computational  problem. 

The  equation  for  the  radial  functions  G  is  obtained  from 
1  n 

Eq.  (4.74)  using  X  (r)=G  fr)/*^.  Equation  (4.74)  is  the  radial 
mm 
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equation  corrected  for  the  presence  of  nonhorizontal  boundaries  and 
is  given  by 


G  ♦  (k  2 (r)  +  -L-  )g  =  -  52  2BmA  *  (5.13) 

m  V  m  4r2  7  m  {feS  mn  n 

The  numerical  procedure  applied  to  the  solution  of  (5.13)  is  the 
iterative  procedure  discussed  in  Chapter  II,  section  C.2,  in  which 
(5.13)  is  replaced  by  the  following  equations  to  be  solved 
iteratively. 


The  zeroth  iterate,  G  ,  is  simply  the  adiabatic  solution. 

m 

Equations  (5.14)  differ  from  (2.57)  and  (2.58)  in  the  omission  of 

the  A  term  which  one  is  not  entitled  to  include  in  first  order 
mn 

calculations,  and  in  the  inclusion  of  a  convergence  acceleration 
parameter,  a,  which  is  discussed  later. 


1.  Numerical  Solution  of  Second  Order  Differential 
Equations 

At  each  step  in  the  iteration  described  by  Eq. 
(5.14)  a  set  of  nonhomogeneous  second  order  differential  equations 
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must  be  solved.  The  numerical  scheme  for  solving  these  differ- 

.  .  .  5. 1-5.4 

ential.  equations  is  Numerov's  method.  *  ’  Numerov's  method  is  a 

multistep,  3-point  computational  scheme  applicable  to  second 

order  differential  equations  of  the  form 

F"  =  UF  +  G  .  (5.15) 

Numerov's  method  is  a  fast  efficient  method  for  solving  second 
order  differential  equations  like  Eq.  (5.15). 

An  equation  of  the  form  (5.15)  is  solved  numerically 
in  Numerov's  method  by  the  following  algorithm. 


(5.16) 


In  Eq.  (5.16)  the  subscripts  refer  to  the  range  mesh  point  at  which 

the  appropriate  quantity  is  evaluated,  and  h  is  the  range  mesh 

increment.  To  get  Numerov's  method  started,  one  must  supply  the 

initial  two  values  of  the  function  F  (F  and  F, ).  With  F  and  F, 

o  1  o  1 

given,  the  solution  can  be  marched  out  using  Eq.  (5.16).  In  the 


150 


calculations  presented  later,  the  initial  values  of  the  radial 

functions  were  obtained  from  the  analytically  known  radial 

functions  for  the  region  r<rg  of  Fig.  V.l. 

The  global^' ^  or  propagation  error  of  Numerov's 
4 

method  is  proportional  to  h  .  The  global  error  is  an  indication  of 
how  much  difference  one  can  expect  between  the  true  solution  of 
(5.15)  and  the  solution  generated  by  (5.16)  at  each  mesh  point.  The 
local  error  of  Numerov's  method  is  proportional  to  h^.  Reference 
5.5  gives  a  detailed  description  of  what  is  meant  by  global  and 
local  error. 

2.  Initial  Conditions 

As  stated  in  the  previous  subsection,  initial  values 
for  the  radial  equations  are  required  to  start  the  Numerov  integra¬ 
tion  scheme.  The  initial  values  of  the  G  functions  of  (5.14)  are 

m 

obtained  by  requiring  the  field  calculated  using  coupled  mode 
theory  in  the  region  r>rg  match  the  field  for  the  range  invariant 
medium  at  the  point  r=rg  (see  Fig.  V.l). 

In  the  source  region,  the  acoustic  field  satisfies 

V2\p  +  =  -  4irqo  6(r-ro)  •  (5.17) 

In  Eq.  (5.17)  the  quantity  qQ  is  the  source  strength  in  units  of 

volume  per  unit  time.  For  the  calculations  to  be  presented  later, 

3 

the  source  was  taken  to  have  unit  source  strength  with  qQ=l  m  /sec. 


The  solution  to  (5.17)  in  the  horizontally  stratified  source  region 
is  given  by  Eq.  (2.21)  which  is 
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<Kr,z)  =  iTTP(z0)  £  Ho1](knr)  VZ)  '  (<i-18) 


In  (5.18)  the  k  and  4>  are  the  eigenvalues  and  mode  functions  for 
n  n 

the  stratified  region  r^T  and  are  independent  of  range  in  that 
region. 

The  initial  values  for  G  used  to  start  the  Numerov 

m 

integration  scheme  are  obtained  by  matching  Eq.  (5.18)  to  Eq.  (5.2) 

at  r  .  If  r  is  the  P„th  mesh  point  in  the  range  grid,  then  G  (r„  ,) 
s  s  n  %- i 

and  G  (r„)  are  taken  from  Eqs.  (5.18)  and  (5.2). 

n  9. 


w 


=  mp 


<v  w  /iyC)(Vj 


(5.19) 


The  6<f>  do  not  enter  (5.19)  since  <5<J>  is  zero  at  r9  ,  and  rp  because 
n  n  "-1  ^ 

the  medium  is  stratified  at  these  two  mesh  points.  The  two  initial 
values  given  by  (5.19)  are  all  that  is  needed  to  start  the  numerical 
Numerov  integration  method. 

In  taking  the  initial  conditions  for  G^  from  (5.18) 
it  is  implicit  that  backscattered  energy  is  being  ignored  as 
discussed  in  Chapter  II.  A  consequence  of  this  assumption  is  that 
the  field  so  obtained  has  only  outgoing  energy  flow.  This  is  the 
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case  because  the  initial  conditions  combined  with  the  numerical 
method  allow  no  mechanism  for  the  generation  of  a  backscattered 
component  of  the  field. 

3.  Convergence  Criteria 

For  the  iterative  process  of  Eq.  (5.14),  some  method 
for  deciding  when  the  process  has  converged  is  required.  There  are 
two  types  of  criteria  that  are  appropriate.  If  Gm  is  less  than  one, 
the  iteration  is  said  to  have  converged  if 


G^1  +  1^  - 
m  m 


<  e. 


(5.20) 


for  all  ranges.  The  factor  e  ^  is  the  absolute  convergence 

criterion.  If  G  is  greater  than  one,  the  iteration  is  said  to  have 
m 


converged  if 


G(i+1)  _  G(i) 

m  m 


•  (i) 


<  e. 


(5.21) 


for  all  ranges.  The  e 2  is  the  relative  convergence  criterion. 
Equations  (5.20)  and  (5.21)  may  be  combined  into  a  single  criterion 
given  by 


cU  +  i)  _  G(i) 

< 

>> 

E.  +  £, 

m  m 

m 

2  1 

(5.22) 
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In  the  calculations  to  be  presented,  Eq.  (5.22)  was 
applied  to  the  real  and  imaginary  parts  of  the  radial  functions  and 
was  required  to  be  satisfied  at  all  range  points  considered. 
Both  and  were  taken  to  be  10  \  With  these  values  for 

and  e2  the  iteration  procedure  is  required  to  converge  to  four 
significant  figures. 

For  the  waveguide  parameters  given  in  Table  II, 
there  are  13  propagating  modes  present  at  the  source  with  mode 
numbers  0  through  12,  and  13  radial  functions  to  be  calculated.  For 
a  bottom  slope  of  5°  a  range  mesh  of  1100  points  was  used.  In  this 
calculation  the  iterative  process  converged  to  four  significant 
figures  after  11  iterations.  The  computer  time  required  for 
convergence  on  the  Applied  Research  Laboratories  Cyber  17]  computer 
system  was  about  367  seconds  or  roughly  6.1  minutes. 

In  Eq.  (5.14)  a  convergence  acceleration  parameter, 
s,  was  included  in  the  iterative  scheme.  In  the  work  presented 
here,  no  convergence  problems  were  encountered  and  s  was  taken  to  be 
zero  in  accordance  with  Eqs.  (2.57)  and  (2.58).  In  some  circum¬ 
stances,  a  nonconverging  iterative  procedure  can  be  made  to 
converge  or  a  converging  procedure  accelerated  by  choosing  s  not 
equal  to  zero  (see  Refs.  5.3  and  5.4).  The  principle  is  similar  to 
the  method  of  successive  over-relaxation  (SOR)  used  in  the  solution 
of  systems  of  linear  equations. 
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C.  Numerical  Results 

In  this  section  numerical  results  obtained  by  the 
corrected  coupled  mode  theory  presented  in  Chapter  IV  are  given. 
The  main  purpose  in  presenting  these  calculations  is  to  illustrate 
the  assertions  made  in  Chapter  IV  and  to  compare  coupled  mode  theory 
with  the  adiabatic  approximation.  Concerning  the  assertions  made 
in  Chapter  IV,  it  will  be  shown  here  that  the  original  coupled  mode 
theory  does  not  conserve  energy  when  applied  to  problems  having 
range  variable  boundaries.  It  will  also  be  shown  that  the  correc¬ 
tions  proposed  in  Chapter  IV  alleviate  this  nonconservation  of 
energy  problem  to  first  order  in  the  bottom  slope. 

1.  Power  Calculations 

In  this  subsection  some  power  flow  calculations  made 
with  the  original  and  the  corrected  coupled  mode  theories  are 
compared.  The  total  power  flow  for  the  calculations  presented  here 
is  taken  in  the  following  form 


where  is  the  power  transported  by  the  nth  mode  in  the  radial 

direction.  In  the  formalism  of  conventional  coupled  mode  theory, 
is  given  by  Eq.  (4.10),  which  is  reproduced  below  in  slightly 
altered  form. 


(n)  _  inco 
r  =  2 


F  F* 
n  n 


'*F  )  +  V  (f 
n  n  /  \  n 


F* 

m 


F*F  )b  l 
n  ml  nm  f 


(5.23) 
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In  obtaining  (5.23)  from  (4.10)  F  =>/r  R  was  used.  The  equation 


n  n 


satisfied  by  F  is 
1  m 


F  +  (k  2  +  -L  )f  =  -  2F  B 
m  \  m  „  2  /  m  n  mn 

4r  n 


(5.24) 


In  the  correction  to  coupled  mode  theory  proposed  in  Chapter  IV  the 
individual  mode  power  amplitudes  are  obtained  from  Eqs.  (4.93)  and 
(4.94)  and  are  given  by 

,S  V 


(n)  =  |/G  g*  -  G*G  )  +  (g  G*  +  G*G  )  ( -r— 

r  2  1  \  n  n  n  n  '  V  n  m  n  m  /  I 

t  V  (g  G*  -  G*G  )  B  l 
\  n  m  n  m  /  nm  [ 
m  j 


mn 


(5.25) 


The  function  G  satisfies  Eq.  (5.13).  In  the  adiabatic  approxi- 
mation  is  given  by  Eq.  (3.16). 

Figure  V.2  shows  a  calculation  of  total  power  flow 
in  the  radial  direction  for  conventional,  corrected,  and  adiabatic 

O  . 

coupled  mode  theory  for  a  bottom  slope  of  5  .  The  source  is  of  unit 

strength  as  defined  by  (5.17)  and  the  power  unit  is  the  watt  (1  W 
2  3 

=  1  kg  m  / sec  ).  The  arrows  along  the  range  axis  denote  the  cutoff 
ranges  of  the  normal  modes.  These  cutoff  ranges  correspond  to  the 
ranges  at  which  (r)  given  by  Eq.  (5.7)  becomes  negative.  Past 
cutoff  the  eigenvalues  become  imaginary  and  the  mode  does  not 
propagate  but  rapidly  decays  exponentially  with  range. 


COUPLED  MODE  THEORY  WITH  CORRECTIONS 
ADIABATIC  MODE  THEORY 


SOURCE  DEPTH  =0,2  km,  c  =  1,5  km/sec 


1  57 


Before  discussing  Fig.  V.2  further,  one  must  know  how  the 

total  power  amplitude  should  behave  in  a  case  in  which  energy  is 

conserved.  It  would  seem  logical  to  expect  the  total  power  flow  to 

be  constant  for  all  ranges  if  energy  is  conserved.  This  is  not  the 

case  depicted  in  Fig.  V.2.  The  reason  for  the  stair-step  type 

behavior  in  Fig.  V.2  is  that  only  forward  going  energy  is  included 

in  the  calculation  of  the  power  amplitudes.  The  stair-step 

behavior  of  the  adiabatic  solution  is  the  expected  type  of  behavior 

for  this  situation.  The  abrupt  drops  in  power  amplitude  occur  at 

the  mode  cutoff  ranges  where  the  outgoing  energy  associated  with 

that  mode  is  suddenly  removed  from  the  problem.  In  between  mode 

cutoff  ranges  the  power  flow  is  constant  indicating  conservation  of 

energy  over  regions  having  a  constant  number  of  modes. 

The  corrected  coupled  mode  theory  calculation  (solid 

curve  in  Fig.  V.2)  which  purports  to  conserve  energy  to  order  H 

should  exhibit  the  same  stair-step  range  dependence.  At  ranges 

greater  than  4  km  this  type  of  behavior  is  present;  however,  at 

lesser  ranges  it  is  not  evident.  This  departure  from  the  expected 

*2  •  3 

behavior  at  short  ranges  is  because  of  terms  of  order  H  and  H  that 

are  introduced  into  the  calculations  through  Eq.  (5.25).  The  G^ 

functions  of  Eq.  (5.25)  satisfy  Eq.  (5.13)  and  contain  terms  of 

order  zero  and  H;  therefore,  since  G^  occurs  quadrat ical ly  in 

•2  '3 

(5.25),  terms  of  order  H  and  H  are  inadvertently  retained.  These 
higher  order  terms  cause  the  departure  from  the  expected  behavior 
at  the  shorter  ranges.  This  departure  is  more  pronounced  at  the 
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shorter  ranges  because  there  are  more  modes  present  at  these  ranges 
and  more  higher  order  terms  introduced  in  the  summation  of  Eq. 
(5.25). 

Notice  the  conventional  coupled  mode  theory  calcula¬ 
tion  (dotted  line)  in  Fig.  V.2.  The  nonconservation  of  energy  is 
immediately  apparent  with  the  observed  rapid  decrease  of  power  flow 
with  range.  Also  notice  the  subtle  differences  between  the 
adiabatic  and  corrected  calculations.  The  effects  of  the  mode 
conversion  process  are  readily  apparent.  The  corrected  coupled 
mode  theory  calculation  shows  an  increased  power  flow  compared  to 
the  adiabatic  calculation  at  the  longer  ranges.  This  increase  is 
due  to  energy  cascaded  from  the  higher  modes  to  the  lower  modes. 
The  power  amplitude  of  mode  zero  (region  to  the  right  of  arrow  1) 
for  the  corrected  coupled  mode  theory  calculation  is  significantly 
larger  than  in  the  adiabatic  case.  The  difference  in  levels  is 
entirely  because  of  energy  transferred  into  mode  zero  from  the 
higher  modes. 

In  Fig.  V.3  the  same  type  of  calculation  as  depicted 
in  Fig.  V.2  for  a  bottom  slope  of  2.5°  is  shown.  With  the  decrease 
in  bottom  slope,  the  higher  order  terms  inadvertently  retained  in 
the  calculation  using  (5.25)  have  a  less  pronounced  effect  as  one 
would  expect.  Once  again,  the  mode  cutoff  ranges  are  denoted  by 
arrows. 

Figures  V.4  through  V.6  are  plots  of  the  individual 
mode  power  amplitude  for  modes  0  through  2.  The  solid  curve  was 


TOTAL  POWER  AMPLITUDE  versos  RANGE 
THE  ARROWS  DENOTE  MODE  CUTOFF  RANGES 
SOURCE  DEPTH  0.2  km,  c  1.5  km/sec 


COUPLED  MODE  THEORY  WITH  CORRECTIONS 


AMPLITUDE  OF  MODE  NUMBER  1  versus  RANGE 

SOURCE  DEPTH  -0.2  km,  c  =  1.5  km/sec 


COUPLED  MODE  THEORY  WITH  CORRECTIONS 
COUPLED  MODE  THEORY  WITHOUT  CORRECTIONS 


SOURCE  DEPTH  =0.2  km,  c  =  1.5  km/sec 
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calculated  from  Eq.  (5.25)  and  the  dotted  curve  from  Eq.  (5.23). 
The  dashed  curve  is  the  adiabatic  power  amplitude  calculated  from 
Eq.  (3.16).  The  adiabatic  calculation  is  missing  from  Fig.  V.6 
because  a  source  depth  of  0.2  km  corresponds  exactly  to  a  null  of  the 
mode  2  depth  function  at  the  source;  hence,  mode  2  is  not  excited 
and  has  no  power  flow  associated  with  it  in  the  adiabatic  approxi¬ 
mation. 

In  Figs.  V.4  through  V.6,  at  ranges  more  than  about 
4  km  there  are  large  differences  between  the  conventional  and 
corrected  theory  calculations.  Once  again,  this  is  a  manifestation 
of  nonconservation  of  energy  caused  by  energy  loss  normal  to  the 
bottom  boundary.  The  effect  is  quite  pronounced,  even  for  a  slope 
of  about  2.5°. 

The  departures  of  the  solid  curves  in  Figs.  V.4-V.6 
from  a  straight  line  are  a  consequence  of  the  mode  conversion 
process.  The  mode  power  amplitude  rises  when  it  is  picking  up 
energy  from  the  other  modes  and  drops  when  it  is  giving  up  energy  to 
the  other  modes.  Figure  V.7  is  a  superposition  of  the  solid  curves 
in  Figs.  V.4-V.6.  Note  the  approximate  negative  correlation 
between  nearest  neighbor  modes.  For  example,  the  peaks  in  the  curve 
for  mode  0  are  roughly  correlated  with  the  valleys  of  the  curve  for 
mode  1.  The  same  type  of  behavior  exists  between  modes  1  and  2  and 
all  of  the  other  neighboring  modes.  This  tends  to  indicate  that  the 
coupling  between  modes  is  strongest  for  nearest  neighbor  modes. 
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Hence,  mode  1  predominately  feeds  and  receives  power  to  and  from 
modes  0  and  2.  This  is  generally  true  for  all  the  other  modes  also. 

2 .  Radial  Functions 

Figures  V.8-V.10  show  the  amplitude  of  the  radial 

functions  G  computed  from  Ea.  (5.13)  for  mode  numbers  0  through  2. 
m 

Also  shown  on  Figs.  V.8-V.9  is  the  adiabatic  approximation  to  |g  I. 
Note  the  differences  between  the  adiabatic  and  coupled  mode 
calculations  that  arise  from  the  mode-mode  interaction  processes. 

In  Fig.  V.10,  no  adiabatic  amplitude  is  included 
because  mode  number  2  is  not  excited  by  the  source;  hence,  in  the 
adiabatic  approximation,  where  no  mode  interaction  occurs,  the 
amplitude  remains  zero  at  all  ranges.  Notice  how  the  mode  coupling 
feeds  amplitude  to  mode  number  2.  Mode  number  2,  which  was 
initially  unexcited  at  the  source,  is  seen  to  build  up  significant 
amplitude  as  it  propagates  away  from  the  source. 

3.  Corrections  to  the  Depth  Functions 

In  this  subsection  the  results  of  the  evaluation  of 
the  corrections  to  the  mode  depth  func  tions  6<}>  ^  are  presented.  The 
differential  equation  for  is  given  by  Eq.  (4.76).  If  the  WKB 

approximation  for  R^/R^  is  employed,  Eq.  (4.76)  becomes 


_£ _  2  2 

a  2  +  k  -  k/(r) 
3z  n 


6<b  =  -  2ik  B  <j> 

yn  n  nn  n 


(5.26) 


The  function  4>  is  given  by  Eq.  (5.9). 
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FIGURE  31.10 

RADIAL  FUNCTIONS  FOR  MODE  NUMBER  2  versus  RANGE 
WITH  1/x/r  RANGE  DEPENDENCE  REMOVED 
SOURCE  DEPTH  0.2  Itm,  c  =  1.5  km/sec 


/ 


As  determined  in  Chapter  IV,  the  desired  solution  of 

(5.26)  is  the  particular  solution.  The  particular  solution  of  Eq. 

(5.26)  may  be  obtained  analytically  in  this  case  from  the  two 
linearly  independent  solutions  of  the  homogeneous  form  of  Eq. 

(5.26) .  The  particular  solution  of  (5.26)  is  given  by  (see  Ref. 

5.7) 


<5<t>n(z;r)  =  -  2ikn(r)Bnn(r)  jy2(z;r)  J' 


z  <t>n(z'  ;r)y1(z'  ;r)dz ' 


A(z') 


-z  <f>  (z'  ;r)y_(z’  ;r)dz' 

-  X^zjr)  J  -B-  2 


A(z-) 


) 


(5.27) 


In  Eq.  (5.27)  y1  and  are  any  two  linearly  independent  solutions 
of 


(5.28) 


and  A(z)  ig  the  Wronskian  given  by 


A(z)  =  y}y^  -  y|y2 


(5.29) 


For  the  case  considered  here  y^  y 2,  and  A  are  given  by 
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y1 (z;r)  =  sin  <n(r)z 
y2(z;r)  =  cos  <n(r)z 


A(z;r)  =  -  k  (r) 


(5.30) 


Kn(r) 


-  kn  ^ 


(n+l/2)u 

H(r) 


Substituting  Eqs.  (5.9)  and  (5.30)  in  Eq.  (5.27)  one  obtains 


60 

n 


nr  k  (r)  B  (r) 
•  n  nn 

1VH  <  (r) 
n 


COS  K  Z 
n 


sin  <  z  3 


(5.31) 


The  solution  for  60  given  in  (5.31)  is  easily  seen  to  satisfy  the 
n 

appropriate  boundary  conditions  given  in  Chapter  IV. 

6<p  (z=0;r)  =  0 
n 


z=ll(r) 


=  ik  H 
n 


.  ( z ;  r) 


z=H(r) 


where  B  is  taken  from  (5.12) 
nn 


B 

nn 


H 
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Figures  V.11-V.13  depict  |  60^  |  versus  depth  for 
mode  numbers  0,  3,  and  5  for  a  range  of  1  km  from  the  source.  The 
bottom  slope  is  2.5  .  Also  shown  on  the  plots  is  the  mode 

depth  function  for  comparison.  Note  that  the  magnitude  of  60^ 
is  indeed  small  compared  to  the  magnitude  of  ^  . 


ODE  FUNCTIONS  <An  AND  8<f>n  FOR  MODE  NUMBER  0 
versus  DEPTH  COMPUTED  FOR  A  RANGE  OF  1  km 
FREQUENCY  =20  Hi,  c  =  1.5  km/sec 


DEPTH  COMPUTED  FOR  A  RANGE  OF  1  km 
FREQUENCY  =  20  Hi,  c  =  1.5  km/sec 


MODE  FUNCTIONS  <£n  AND  S<£n  FOR  MODE  NUMBER  5 
versos  DEPTH  COMPUTED  FOR  A  RANGE  OF  1  km 
FREQUENCY  =20  Hz,  c  =  1.5  km/sec 
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At  this  point  it  should  be  emphasized  that  6$  may 
be  determined  analytically  in  this  case,  only  because  $  is  known 
analytically.  In  a  more  realistic  problem  in  which  the  must  be 
computed  numerically,  the  6$^  will  also  have  to  be  computed 
numerically. 


4.  Field  Calculations 

In  this  final  subsection  of  Chapter  V  results  of  the 
calculation  of  the  acoustic  field  iJj  are  given.  When  considering 
the  acoustic  field,  it  is  customary  to  examine  \p  as  a  function  of 
range  for  specific  source  and  receiver  depths.  The  calculations  to 
be  presented  here  were  made  for  the  waveguide  model  of  Fig.  V.l  with 
the  source  located  at  a  depth  of  0.2  km  and  the  receiver  at  a  depth 
of  0.1  km.  The  slope  angle  was  taken  to  be  5  . 

Rather  than  examining  4>  as  a  function  of  range,  the 
absolute  square  of  *P,  is  considered,  as  is  customary  in  under¬ 

water  acoustic  modeling  work.  The  absolute  square  of  the  field  is  a 
more  meaningful  quantity  since  it  is  proportional  to  the  measured 
intensity  given  by 

P2 

rms 

I  =  -  •  ( s  T9 1 

meas  p  c 
o  o 

In  Eq.  (5.32),  p  is  the  root  mean  square  pressure  measured  at  the 
rms 

receiver,  and  Pq  and  cq  are  the  density  and  sound  speed  at  the 

receiver.  The  constant  of  proportionality  linking  *meas 

and  is  w2  p  /2c  . 

o  o 


(5.33) 
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Figure  V.14  shows  ip*tp  versus  range  for  the  indicated  source- 
receiver  depths.  Two  curves  are  presented  in  Fig.  V.14,  one  for  the 
adiabatic  approximation  of  and  one  for  a  coupled  mode  theory 

evaluation.  Notice  how  the  fine  structure  of  <p*  differs 
appreciably  between  the  adiabatic  and  coupled  mode  calculations. 

Figure  V.15  shows  the  difference  between  the  two 
curves  of  Fig.  V.14  as  a  function  of  range.  It  is  obvious  from  Fig. 
V.15  that  some  large  differences  exist  between  the  adiabatic  and 
the  coupled  mode  calculations.  If  the  curve  in  Fig.  V.15  is 
averaged  over  range,  one  finds  the  mean  to  be  +0.083  which  is  very 
close  to  zero.  Hence,  the  range  averaged  difference  between  the 
adiabatic  and  coupled  mode  calculations  is  very  slight.  In  many 
underwater  acoustics  applications,  the  fine  structure  of  the  field 
is  not  of  interest.  Many  data  analysis  procedures  involve  a  range 
averaging  operation  on  the  measured  field.  It  seems  plausible 
that,  in  situations  such  as  these,  an  adiabatic  calculation  of  the 
field  might  provide  an  adequate  characterization  of  the 
medium,  which  would  greatly  simplify  the  computational  requirements. 

Figure  V.16  shows  the  phase  structure  of  the  field 
used  to  generate  Fig.  V.14.  The  adiabatic  and  coupled  mode  theory 
calculations  of  the  phase  are  both  depicted.  The  reason  for 
including  this  figure  is  to  check  the  assertion  made  in  Chapter  III 
concerning  the  phase  characteristics  of  the  adiabatic  solution.  In 
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FIGURE  31.15 

DIFFERENCE  BETWEEN  COMPUTED  BY  COUPLED  MODE  THEORY 
AND  ADIABATIC  MODE  THEORY  versus  RANGE 
SOURCE  DEPTH  =  0.2  km,  RECEIVER  DEPTH  =  0.1  km.  c  =  1.5km/*.e 
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FIGURE  2.16 

PHASE  OF  FIELD,  *,  IN  RADIANS  versus  RANGE 
SOURCE  DEPTH  =0.2  km,  RECEIVER  DEPTH  =0.1  km,  c  =  1.5  kmAec 
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Chapter  III  it  was  concluded  that  the  adiabatic  phase  of  the 
acoustic  field  should  be  a  good  approximation  to  the  phase  of  the 
field  computed  using  coupled  mode  theory  with  the  mode-mode 
coupling  effects  primarily  influencing  the  scalar  distribution  of 
amplitude  among  the  modes.  Figure  V.16  supports  this  claim.  The 
adiabatic  phase  tracks  the  coupled  mode  phase  quite  well.  In  fact, 
the  agreement  is  better  than  indicated  by  Fig.  V.16  since  a  good 
portion  of  the  difference  between  the  two  curves  is  due  to  sampling 
effects  and  the  resolution  of  the  plotting  hardware  used  to 
generate  Fig.  V.16.  The  finding  that  the  adiabatic  phase  is  a  good 
approximation  to  the  true  phase  of  the  field  is  important  because  in 
many  applications,  such  as  modeling  array  response,  the  phase  of 
the  field  is  of  more  importance  than  the  amplitude  information. 
Hence,  it  seems  likely  that,  in  cases  such  as  this,  an  adiabatic 
characterization  of  the  field  would  be  adequate,  again  greatly 
simplifying  the  required  computational  effort. 


CHAPTER  VI 


CONCLUSION 

The  purpose  of  this  report  has  been  to  examine  the 
theory  of  sound  propagation  via  coupled  normal  modes  as  applied  to 
acoustic  propagation  in  a  range  dependent  ocean  environment.  The 
emphasis  has  been  placed  on  low  frequency  propagation  over  range 
variable  ocean  bottoms.  Coupled  mode  theory  and  the  adiabatic 
approximation  to  coupled  mode  theory  involve  approximations 
concerning  the  boundary  conditions  satisfied  by  the  field  and  the 
strength  of  the  range  dependence  of  the  medium  respectively. 
Before  the  theory  can  be  used  as  an  acoustic  propagation  modeling 
tool,  the  effect  of  the  approximations  involved  in  the  theory  must 
be  examined  with  regard  to  the  types  of  range  variability  expected 
to  be  encountered  in  realistic  ocean  media.  This  examination  was 
the  main  purpose  of  the  work  presented  in  this  report. 

In  Chapter  II  the  general  framework  of  the  theory,  as 
originally  proposed,  was  described.  In  the  derivations  all 
quantities  which  must  be  evaluated  in  the  implementation  of  the 
theory  were  given  in  a  form  suitable  for  numerical  computation.  The 
methodology  for  implementing  the  theory  in  realistic  problems  was 
also  described. 

Chapter  III  considered  the  adiabatic  approximation  to 
coupled  mode  theory.  Two  aspects  of  the  adiabatic  approximation 
were  investigated:  (1)  its  sensitivity,  as  a  function  of  sediment 
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type,  to  radial  sound  speed  derivative  and  water-sediment  interface 
slope,  and  (2)  whether  the  geometric  (raypath)  properties  of  the 
adiabatic  field  were  given  correctly.  The  adiabatic  approximation 
was  seen  to  be  sensitive  to  bottom  type  with  respect  to  lateral 
sound  speed  gradient  and  water-sediment  interface  slope,  and  the 
nature  of  the  sensitivity  was  explained  using  ray  theory  arguments. 
The  adiabatic  solution  was  also  seen  to  describe  the  raypath  and 
phase  structure  of  the  field,  even  for  quite  large  rates  of  radial 
change.  This  finding  led  to  the  conclusion  that  the  mode-mode 
coupling  effects  neglected  in  the  adiabatic  approximation  were 
primarily  responsible  for  amplitude  redistribution  leaving  the 
phase  structure  essentially  unchanged. 

Chapter  IV  examined  the  boundary  condition  approximation 
inherent  in  coupled  mode  theory  when  applied  to  problems  with  range 
variable  boundaries.  Since  nonhorizontal  boundaries  are  present  in 
almost  all  ocean  propagation  problems  of  interest,  the  work 
presented  in  this  chapter  is  probably  the  most  significant  portion 
of  the  dissertation.  In  Chapter  IV  the  consequence  of  the  boundary 
condition  approximation  was  explored.  It  was  found  that  the  order 
of  the  approximation  on  the  boundary  conditions  was  not  consistent 
with  the  order  of  coupled  mode  theory.  It  was  also  shown  that  a 
consequence  of  this  inconsistency  was  energy  not  being  conserved 
within  the  original  formulation  of  coupled  mode  theory.  A  con¬ 
sistent  correction  to  coupled  mode  theory  was  then  derived  such 
that  the  field  so  obtained  satisfied  the  proper  boundary  conditions 
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and  conserved  energy  to  first  order  in  the  local  slopes  of  any 
nonhorizontal  boundaries.  The  corrections  to  the  theory  were  given 
in  a  form  convenient  for  computational  purposes.  The  results  of 
Chapter  IV  have  exposed  the  consequences  of  the  inconsistent 
boundary  condition  approximation  inherent  in  coupled  mode  theory, 
heretofore  ignored,  and  have  shown  how  to  overcome  the  inconsis¬ 
tency. 

Chapter  V  presents  the  results  of  numerical  calculations 
that  were  performed  for  the  purpose  of  illustrating  and  verifying 
the  assertions  and  corrections  proposed  in  Chapter  IV.  The 
numerical  calculations  given  in  this  chapter  were  made  for  an 
isovelocity,  wedge  shaped  waveguide  having  perfectly  reflecting 
boundaries.  Power  flow  calculations  were  presented  to  illustrate 
nonconservation  of  energy  in  the  original  formulation  of  the  theory 
and  to  show  that  the  proposed  corrections  produced  the  desired 
results.  Calculations  of  the  radial  and  depth  functions  as  well  as 
a  calculation  of  the  field  were  presented. 

There  are  several  possible  extensions  of  the  worn 
presented  in  this  report.  One  concerns  the  inclusion  of  solid 
properties  of  the  medium.  Coupled  mode  theory  is  based  on  the 
assumption  of  a  fluid  acoustic  medium.  In  realistic  oceanic  wave¬ 
guides  the  fluid  assumption  is  valid  for  the  water  layer  and  is  a 
good  approximation  for  the  unconsolidated  sediments.  However,  the 
underlying  basalt  substrate  which  occurs  at  varying  depths  beneath 
the  ocean  bottom  is  known  to  exhibit  significant  solid  properties. 
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Hence,  one  possible  extension  of  the  work  given  here  would  be  a 
reworking  of  the  theory  to  include  the  solid  properties  of  the 
basalt.  Such  an  extension,  however,  would  likely  make  an  already 
complicated  formalism  prohibitively  complex. 

Another  possible  extension  concerns  the  importance  of  the 
nonconservation  of  energy  effect  compared  with  any  realistic 
attenuation  processes  occurring  in  the  medium.  This  comparison  was 
suggested  in  Chapter  IV.  In  the  numerical  calculations  presented 
in  Chapter  V  the  effect  of  nonconservation  of  energy  was  seen  to  be 
severe.  It  might  turn  out  that,  in  more  realistic  waveguides  with 
attenuation  included,  the  effect  is  negligible.  If  this  were  true 
in  some  cases  the  implementation  of  the  theory  could  be  con¬ 
siderably  simplified  because  the  corrections  proposed  in  Chapter  IV 
could  be  ignored. 

A  third  possible  extension  concerns  the  finding  of 
Chapter  III  that  the  true  phase  structure  of  the  field  is 
essentially  given  by  the  adiabatic  approximation  to  coupled  mode 
theory.  A  consequence  of  this  finding  is  that  in  applications  of 
the  theory  in  which  the  phase  structure  of  the  field  is  the  primary 
concern,  e.g.,  array  processing  and  field  coherence  studies,  the 
adiabatic  approximation  is  perhaps  all  that  is  required  to  obtain 
meaningful  results.  Since  this  would  involve  a  great  simplifi¬ 
cation  of  the  calculations,  it  is  worth  examining  further. 


An  obvious  extension  of  the  work  presented  in  Chapter  IV 
is  to  include  second  order  corrections  to  the  field  and  boundary 


conditions  and  estimate,  in  terms  of  bottom  slope,  when  their 
inclusion  is  necessary.  Such  an  extension,  however,  would  probably 
not  be  computationally  feasible  because  of  the  added  complexity  it 
would  entail. 
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SOME  PROPERTIES  OF  THE  NORMAL  MODE  EIGENFUNCTION  EXPANSION 


The  purpose  of  this  appendix  is  to  present  a  more 
detailed  description  of  the  properties  of  the  normal  mode  depth 
function  described  by  Eq.  (2.12)  which  is  reproduced  below. 
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The  boundary  conditions  on  4>n  are 
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The  sound  speed  structure  of  a  typical  oceanic  waveguide  is  such 

that  the  wave  number,  k(z),  attains  some  limiting  value  as  z  goes  to 

infinity,  k  .  .  The  maximum  value  of  k(z)  occurs  in  the  water 
min 

column  at  the  sound  speed  minimum  which  defines  the  center  of  the 

sound  channel.  Call  this  value  of  k,  k  .  Equation  (A.l)  may  now 

’  max  ^  J 


be  rewritten  as  follows. 
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4>n(z)  =  0 


(A. 3) 
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The  function  q(z)  as  defined  by  (A. 4)  is  negative  and 
appioaches  zero  as  z  goes  to  infinity.  Equation  (A. 3)  and  boundary 
conditions  define  an  eigenvalue  problem  whose  properties  have  been 
widely  investigated  (see  Refs.  A.l  and  A. 2).  The  eigenvalue 
spectrum  of  Eq.  (A. 3)  is  known  to  be  discrete  for^^O  and  con¬ 
tinuous  for  ^n>0.  For  readers  familiar  with  quantum  mechanics,  the 
problem  defined  by  Eqs.  (A.3)-(A.5)  is  analogous  to  the  Schrodinger 
equation  for  a  particle  in  a  potential  well  where  particles  with 
energies  less  than  zero  occupy  bound  states  with  discrete  energy 
eigenvalues  and  particles  with  positive  energies  can  have  a 
continuum  of  energy  eigenvalues.  Therefore,  in  the  notation  of  Eq. 
(A.l)  the  eigenvalues  are  discretely  distributed  for  magnitudes 


— ajgrgf ■  m 


187 


between  k  .  and  k  and  are  continuously  distributed  between  0 
min  max 

and  k  .  . 
min 

The  mode  functions  $  including  the  continuum  form  a 
complete,  closed  set  of  functions  on  the  depth  interval  (0,°°) 
having  the  following  closure  property. 


n=l 


p(z)<HK,z)<KK>z')dK  =  6(z-z')(A.6) 


In  (A. 6)  the  4>(<,z)  is  a  continuum  mode  of  eigenvalue  K.  The 
equation  for  the  field  from  a  point  source  given  in  Eq.  (2.21)  with 
the  inclusion  of  the  continuum  is 
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The  depth  functions  satisfying  Eq.  (A.l)  are  orthogonal 
with  respect  to  the  density  as  a  weighting  function.  To  see  this 
consider  the  following  expression  obtained  from  Eq.  (A.l). 
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If  Eq.  (A. 8)  is  multiplied  by  the  density  which  is  piecewise 

continuous  throughout  the  layers,  and  the  result  integrated  over  z 

the  following  is  obtained, 
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(kn2  -  \2)f  p(z)<t>„(z)4>Jz)dz  =  0 


m  n 


(A. 9) 


The  RHS  of  (A. 9)  is  zero  because  of  the  boundary  conditions  given  in 
Eq.  (A. 2).  Equation  (A. 9)  shows  that  the  mode  functions  of 
different  mode  index  are  orthogonal  to  each  other.  The  same  argu¬ 
ments  show  that  discrete  and  continuous  modes  are  orthogonal  to 
each  other  as  are  two  continuous  modes  with  different  continuum 
eigenvalues.  The  density  p(z)  plays  the  role  of  a  weighting 
function  in  the  orthogonality  relationship.  It  is  convenient  to 

normalize  the  <f)  such  that 
n 


f  P(z)<f>n(z)<j>m(z)dz  =  Sn>) 


(A. 10) 


The  continuous  modes  are  normalized  as  follows: 


f 


p(z)<Hk,  z)4>Ck'  » z)dz  =  5(k-k') 


(A. 11) 


In  Chapter  II  it  was  stated  that  the  contribution  to 
the  field  from  the  continuum  is  negligible  in  many  cases.  To  show 
that  this  is  true,  a  proof  given  by  Mitchell  is  reproduced. 
Consider  the  continuum  contribution  to  the  field  in  Eq.  (A. 7). 


tl>  =  iTTpfZ  ) 

rcont  oJ  j i  o 
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(A. 12) 


Divide  the  interval  (0,k  .  )  in  such  a  way  that 

mm 


£  ¥(*) 


(A. 13) 


ip(SL)  =  iirp(z 


(Kr)c|>(K,z)<t>(K,  z0)dK  .  (A.  14) 


The  intervals  defining  the  ’P(^)  are  chosen  such  that  IT  (<r) 

o 

vanishes  at  both  end  points  and  such  that  <{>  (k ,  z)4>  (k  ,  z^)  has  the 
same  sign  over  the  interval.  With  this  choice  of  intervals  one  can 
set  an  upper  bound  on  the  M&), 


-  i  r*+i  a) 

^  ttG CIO  I  HLU 

1  J  „  ° 


(Kr)  d< 


I  *1 


(A. 15) 


*(*)  *^1/ 


(x)dx 


(A. 16) 


G(H)  =  max[4>(K,z)<)>(ic,zo)] ;  <  tc  <  k£+1 


If  the  asymptotic  form  of  the  Hankel  function,  which  is  valid  for 
ranges  more  than  a  few  wavelengths  from  the  source,  is  used  in  Eq. 
(A. 16),  one  obtains 


<  const 


i 


Tcont 


<  const 
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The  field  from  the  continuous  modes  is  therefore  seen  to  fall  off  as 
l/r  '  for  ranges  more  than  a  few  wavelengths  from  the  source.  This 
falloff  is  l/r  times  faster  than  the  l//r  range  dependence  of  the 
discrete  contribution  to  the  field.  Hence,  one  can  see  that  the 
continuous  mode  contribution  to  the  acoustic  field  is  negligible  in 
the  farfield.  At  ranges  in  close  proximity  to  the  source  and  in 
problems  involving  only  a  few  discrete  modes,  the  continuum  can  be 
important.  Papers  by  Stickler, Williams,^"  and  Tindle  *  have 
examined  these  situations. 
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